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I. INTRODUCTION 


The importance of vortical flows in a variety of aerodynamic applications i3 well-known 
and has been the subject of study for almost a hundred years. Despite numerous theoretical 
and experimental studies of these flows, however, they are still not fully understood and 
therefore still continue to gain the attention of the aerodynamicist. In this report three 
selected problems involving vortical flows are analysed and discussed, namely 

(i) The trailing vortex behind a wing 

(ii) Rotor blade-vortex interaction 

(iii) The leading edge vortex on a flat plate 

On each case a simplified model of the flow is formulated and approximate solutions 
are found to gain insight into the basic physical phenomena involved. 

n. DISCUSSION OF SELECTED PROBLEMS 
(i) Trailing Vortex 

Vortices appear as a necessary consequence of aerodynamic lift and persist indefinitely, 
in the absence of laminar or turbulent diffusion, after the passage of an aircraft. In this 
regard their rate of decay due to diffusion is of practical interest inasmuch a3 their continued 
presence behind an aircraft constitutes a potential hazard to following aircraft particularly 
at low altitude in the terminal area. An analysis of the decay of a vortex pair generated 
by a wing is given in Appendix A. 

The analysis provides approximate closed form solutions for each of two regions (a) a 
persistence region in which the circumferential flow within the vortex remains unchanged 
but the axial flow accelerates to the free stream value, and (b) a decay region in which the 
circumferential flow decays through turbulent diffusion while the axial flow remains at the 
free stream value. 

In particular the analysis provides expressions for the length of persistence of the 



vortex before the onset of decay in terms of the aircraft span, wing aspect ratio and lift 
coefficient. The analysis also gives the vortex radius and velocity profile as a function of 
distance behind the aircraft. 

(ii) Rotor Blade- Vorter Interaction 

One of the primary sources of helicopter noise, under certain flight conditions, arises 
from the interaction of the vortex produced by a blade with the aerodynamics of the 
following blades. This phenomena, sometimes called blade-slap, 'is discussed and analyzed 
in Appendix B. The analysis has three major parts (a) the rotor-tip vortex, (b) the unsteady 
aerodynamics of the interaction with the following blade, and (c) the acoustic pressure field. 

The analysis of the vortex formed by the blade tip includes the effects of turbulence and 
viscosity at the center of the vortex and therefore does net introduce a singularity as is the 
case for inviscid vortex models. In the unsteady aerodynamics analysis the disturbance due 
to the passage of the vortex over the blade is represented as an unsteady gust using a Kutta 
condition at the trailing edge and the assumption of a flat wake. Acoustic pressures arc 
then determined, using the Hawkings-Ffowcs Williams equation and the results expressed 
in terms of the vortex characteristics, ie vortex strength, core size and Reynolds number. 

A comparison of vortex velocity profiles are found to be in good agreement with exper- 
iment. The unsteady lift is found to be dependent on the core size (which in turn varies 
as the blade radius), on distribution of blade loading and on the Reynolds number based 
on vortex strength. Extensions of the analysis to determine the effect of spanvvise uon- 
compact sources and the effect of variations in blade tip chape on acoustic signature were 
also undertaken. There appears to be a significant difference between the results of the 
non-compact source analysis and those for the compact source; even in the far fields there 
is an appreciable variation in acoustic signature, both in amplitude and shape, depending 
on observer position. Finally the effect of tip chape on acoustic signature is determined 
and compared to the effect of tip shape on blade aerodynamics efficiency. Reductions in 



acoustic pressure of as much as 6dB may be realised by tip shaping but with an appreciable 
loss of lifting efficiency (of the order of 20% - 30%). 

(iii) Leading Edge Vortex on a Flat Plate 

The classical solution for flow over a fiat plate in 2-dimensions satisfies a trailing edge 
Kutta condition but produces a singular behavior at the leading edge. An alternative 
model has been studied which requires the introduction of a vortex sink placed above the 
leading edge of the plate. The flow equations are solved and the position and strength 
of the vortex sink determined in terms of the angle-of-attack of the plate. The eolation 
so found represents the flow near the root of a thin wing where the sink gives rise to a 
strong rpanwisc flow contained within the detached vortex, and gives a non-zerc drag for 
the plate. 

A comparison is made of this solution with the results of two classical solutions, namely 
the attached flow with the leading edge singularity and the totally separated flow (Helmhob 
Solution). The vortex- sink solution is found to exist in the angle of attack range 0 < a < 
72° with maximum lift occurring at a = 45*. Thi3 work is given in detail in Appendix C. 

m. CONCLUDING REMARKS 

The understanding of vortical flows continues to need further improvement in view of 
the fundamental role they play in a variety of aerodynamic problems and applications. 
Although much progress has been made in the analysis of inviscid vortical flows there con- 
tinues to be a lack of complete understanding in areas involving vortical flows interacting 
with turbulent phenomena. In some cases, as for example in the case of the trailing vortex, 
such interactions can not be ignored since they determine the rate of vortex decay (and 
therefore the extent of the hasard to a following aircraft). In other cases, as in the case 
of the blade vortex interaction, turbulence removes what would otherwise be an inviscid 
singularity in the pressure at the vortex center, and permits a realistic assessment of the 



acoustic signal propagated by the blade. Other problems such as vortex bursting, which 
lie outside the scope of this report are also of interest in aerodynamic design and deserve 
increased attention in the future. 
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SUMMARY 

This paper provides an analysis that describes the trailing vortex generated by a lifting 
surface, the structure of its turbulent core and the influence of axial Bow within the vertex on 
its initial persistence and on its subsequent decay. Similarity solutions of the turbulent diffusion 
equation are given in closed form and results are expressed in sufficiently simple terms that the 
influence of the lifting surface parameters on the length of persistence and the rate of decay of 
the vortex can be evaluated readily. 
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I. INTRODUCTION 


SYMBOLS 
wing reference area 
aspect ratio b’/A 
aircraft wing span 
distance between vortices 
induced drag coefficient 
lift coefficient 
2rr/* : r, 

characteristic length of vortex persistence 
wing lifting efficiency. C\/zAHCp 

= h 

constant. 0 06 

outer edge of turbulent core 
L'b/v.Ul 

distance in radial direction 

lift distribution parameter. J 0 f-d)‘ 

=sinh(4s : /e - 11/12) 

velocity function 
maximum value of V 
circumferential velocity 
distance in axial direction 
r/d 

=y/(6/2) 

distance in spanwise direction 
— r /r,(x) 

circulation 

coefficient of eddy viscosity 
coefficient of kinetic viscosity- 
density 


The structure of the wake vortex generated by a lifting surface «ucb as an aircraft wing or 
a rotor blade has been the subject cf much previous study, both theoretical and experimental 
(see for example ref. 1-t). The trailing vortex, after its initial formation is known to comprise 
two phases, namely (!) a region of persistence for long distances compared to the wing spaa and 
(2| a region of decay in which the circumferential velocity is reduced, and the vortex diffuses 
radially with increasing distance (fig. I|. In this regard it is of particular interest to understand 
the structure of the vortex in sufficient detail to predict the velocity distribution in the vortex, 
to explain the initial persistence, and characterize its subsequent rate of decay in terms of the 
aerodynamic properties of the lifting surface 


| 
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2. ANALYSIS 


The I\>tential Flow 


The potential flow due to a vortex pair his been analyzed by Spreite; and Sachs (ref l). In 
their approach the vorticea are considered as a pair of rotating cylinders around which a potential 
flow exists in the plane normal to the axis of the vortices. The relevant expressions which relate 
the characteristics of the wing to those of the vortex pai. are derived here, for completeness, as 

follows. 


The lift is expressed as 


L 


U 2 C L A 



(II 


where T| is the circulation at the root, G «=* is the spanwise distribution of lift, and Y «= 
is the dimensionless spanwise distance. 


The integral in equation (1) is denoted by 


£ 


GiY ~ i 


j is a parameter which characterizes the lift distribution. 


( 2 ) 


The vortex sheet behind the wing rolls into a cylindrical vortex, the inner part of the cylinder 
containing the vorticity shed from the tip and the outer part containing that shed from the root; 
thus T| is the strength of the vortex at its outer edge r t . 

From equation (1), T| may be written 


r ~r b £L 1 

P| L 2 Alt' s 


where the aspect ratio is Ap = b'/A. 


(3) 


Following further the Spr»iter-Sachs analysis, the induced drag is Mated to the rate of 
formation of rotational kinetic energy in the wake. This relationship is written 

D = I^CdA = ff '-px-idS' (4) 

where the integral is performed over the entire plane normal to the direction of flow. 


Equations (I) and (1) can be combined to given an expression which relates the energy integral 
to the wing efficiency 


//(f0«'= T 


( 5 ) 


where e = rtjfcp is the wing efficiency factor. 


The integral in equation 3 has contributions from the potential flow outside the vortex pair 
and from the vortex pair itself; thus 


//(sM//teH 4m'T 

JfOl 


( 6 ) 


where the substitution r* = his been made within the vortex pair. 


The contribution from the potential Row is evaluated in closed form and may be written 



T 


3-3 


uhtrf t] is the separation distance between the vortices. The contribution from the vortex pair 
must be evaluated from the radial distribution cf £ = C within the vertex. This contribution 
is found by solving the turbuhnt Cow equations and. as seen later, gives 



Substitution of (6). (7), and (8) into (5) gives an expression for n, 



(8) 


(0) 


where S' == ?inh( ts : /e — 11/12) 


The relationship (0) between the vortex radius r-| and the separation distance t| will hold 
even when n varies with the distance z along the vortex since the flux of rotational energy must 
be constant and equal to the induced drag at all values of z. 


A determination of the separation distance t| between the vortices requires that an additional 
condition be prescribed. It is assumed that the centroid of vorticily is conserved during the roll-up 
process and thus the initial spanwise separation of the vortex ts given by 


MOI 



(after integration by parts). 

The initial radius of the vortex is 


ri(0) = (10) 

The initial vortex radius rj(0) and the circulation T|, given by equation (3) are important 
parameters in the description of the turbulent core of the vortex ej seen in the following analysis. 


The Turbulent Yor'ex Core 


For any location z downstream of the lifting surface the characteristics, of the vortex will be 
influenced by turbulence and. to a lesser degree, by kinematic viscosity. There will exist a core 
within which the circumferential velocity will differ from that given by the potential flow (fig 
21. This vortex core flow can be described by a differential equation expressing the transport of 
angular momentum as follows: 


fG tG 
«r + v— 
6z er 


. 

t- + < - — 

r 6r 


(111 


where v is the kinematic viscosity and i is a turbulent eddy viscosity. 


The eddy viscosity for vortex flows has been investigated by Hoffman and Joubert (ref 5) and 
a model formulated in a way analogous to the traditional mixing length theory of I’randtl. A 
dimensional analysis suggests that the eddy viscosity is related to the circumferential shear stress 
through the relationship 



r 



which may be comitined to give 


where the constant k must he determined experimentally. The value k = 0 03 appears to fit the 
experimental data. 

0 




3-4 


When equation (12) is substituted into equation (11) tht ."Hull is written 


(C IG k'-T, S \f fG\l(G 
* 6 1 + * 6r = 2i r tr\+ r 6r)r6r 


(13) 


2tv irv.iR I 

f = W, “ FTTcjT 

(substituting for T| from equation (3)1. 

The quantity e is seen to vary inversely with the Reynolds number based on a chatac- 
teristic chord length. -£p. Tor a typical aircraft wing e is of the order of 10' 2 or 10~ 5 and the 
efTec's of kinematic viscosity on the overall wake characteristics are negligible. However, la tee 
region close to the center of the vo’tex (where r £ .G < r) tht effect of viscosity will change the 
local velocity gradients significantly, creating a "laminai sub-core*. 

Equation 13 must be solved for the distribution (7(r) subject to appropriate boundary condi- 
tions, and using suitable approximate expressions for u and f. the convective velocity components 
Two cises are considered here corresponding to the region of persistence, where the vortex remains 
tightly rolled, and the region of decay where the vortex increases its radius and decreases its 
circumferential velocity with distance along the vortex. 


(I| Region of Persistence : Immediately behind the wing it is to be expected that the axial velocity 
u in the vortex co'e will be less than the free stream value, V but must accelerate to this vah.e in 
a distance d «av. The velocity components along and uormal to the centerline of the vortex may 
be approximated as 

u I 

r = d' 

satisfying the equation of continuity. 


* 1 T , 

r =--- Jot red 


The accelerating flow along the axis (see fig. 3a) causes a radial inflow which in turn converts 
vorticity inward to balance the outward diffusion by turbulence. This causes the vortex to remain 
tightly rolled until the axial velocity returns to its free-streara value and the inflow ceases. 

(Jl Reginn of Decav : A*, far distance behind the wing where the axial velocity has recovered to 

the free siream value 

— s= l . — ssv 0 for z > d. 

Here (see fig 3b) there is no radial inflow but radial diffusion is balanced by axial convection 
resulting in » spreading of the vortex and as will be seen, a decay in the circumferential velocity, 
with distince along the vortex. 

Approximate Solutions 

For both of these regions simple approximate solutiona for G can be obtained from equation 
(II). First it is convenient to define dimensionless variables 


7 t u ri(x) _ 

l = - : -- ; - ss l*. -id-J = »"> 

d r,(x) L fi(l) 

and seek, self similir solutions of the form G «= f7(r) 

The resulting ordinary differential equation, derived from equation 13, is written 

t 

+ iJ-G 1 = o 


[(<• + sC>— 


where 


and the prime denotes differentislion with respect to g. 


[ gf’rj(l) 

d 


: um] 

L * : r,d 

di 

x r i i 1 1 / 

‘■J 


(14) 


(IS) 


Sin"e ,f must be independent of 1 for similarity to hold, the exponents m and n in the 
express)' ns for r | ( J ) and u(l) respectively must satisfy the relation 2m + n = 1. thus giving from 
equation (15). 

J 5 =3 irfT'(l)/Jt 2 r,d (15a) 
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ORIGINAL PA C;Z 
OF POOR QUALITY 


Since m = it can be seen that ia the region of persistence ( J < 1. ^ * n = l) the 

corresponding value of m is »n = 0, so that ri(z) - r|(l) = ri|0). Similarly, in the region of 
decay (J > l,^r = l,n=0)tne value of m is m = I, so that r^r) = n ( 1 )( *)■ . 

In order to obtain the solution. C(z). equation (11) must be solved subject to the boundary 
conditions 1 > = l.C'(l| = 0 (assuming that the core merges with thp potential solution (G = 
II. ) and GlO) = 0 giving zero circumferential velocity at the center of the vortex. 


Equation 14 is nonlinear; however an approximate solution can be found in two steps: since 
t is a small parameter ihe solution G for c = 0 is found valid except for a small region near 
the center of the vortex where zG' < c, and an improved solution, using the approximation 
c + :G' = c + zG , is then found from the resulting lineat equation. 

Thus, with c = 0, equation 14 becomes 



+ 4 3 i (/ = 0 


(14a) 


and a solution satisfying the boundary conditions is determined as 


I- G = (l---) 5 
3 = 1 

The term (c + zG') is now approximated as 

r + zG* « c+ zG' s; 2(1 - z)(z + £) 


( 10 ) 


so that equation 14, with 3 = 1, becomes 


f O' 


+ 2 (f - 0 


Which has the improved solution 


1 -C = <1 
3 = 1 



( 116 ) 


(IGa) 


Equation IG closely approximates equation (16a) for small values of c except for z = 0(c), 
corresponding to a laminar subcore' at the center of the vortex. 


All of the relevant characteristics of the vortex can now be determined using equation (16a), 
i.e., (he length of persistence, the radius of the turbulent core, the radius of the laminar sub-core, 
the velocity profiles and the variation of vortex size and velocity with distance i along the axis. 


First the integral of equation 8 can now be evaluated using the approximate solution G = 
I - (1 - z| 3 giving 



111 
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thus verifying equation 8, 

The characteristic distance of persistence , d is found from 15a, with 3 — 1, as 


d 


* V 

T k-r i ’ 


( 17 ) 


The same result, expressed in terms of the aerodynamic parameters of the wing {6.C/., and AR) 
is written 




(17a) 


The radius of the turbulent core rj is determined as 


r ° r T < 1 
a 


n = t-i(o) 

...rrt! 


( 18 ) 



JbI2&A£tl ?^£J.±2-l w.'x. - -.. 


corresponding to the region of persistence and decay respectively. 


Alternatively, in terms of the wing parameters 


n = «$- 


for | < 1 


fSk 1 C t \* _i fz\\ b 1 

“VT-^sJ* *Uh r ° r d>‘ 


with d given by equation 17a. 

The circumferential velocity is determined from <> = as 
*’* == 2^) V ' , ** K *d <l 


= s£ip» v Kf)" t - 3 > 1 


a i-u-*) s (i+¥r 

'•(-•) = 7 = ; . * < i 

*->i 


with a = and ri(r) given by equation (18), d by equation (17). 


Alternatively 




for - < 1 

a 


,oti d >i 

with d given by equation (17a). The variation of v$ with j is shown in fig. 4. 

The velocity profile V(a) represented by equation 20 has the following characteristics (described 
in fig. 5): 

(a) For 2 <C c, V'(r) a; "f giving a linear variation with slope * near the center of the 
vortex corresponding to the laminar subcore, 

(b) For c <3C 2 < 1; l'(r) = 2 — t corresponding to the turbulent core, 

(c) For 2 > l;\’(r) = J corresponding to the potential vortex. 

The maximum circumferential velocity is found by differentiation of V(r). For c < 1 it can 
be shown that the miximum occurs at £ = a* = (S log and has the value 


f-s -agio *!) 1 


Tor most purposes, the limiting case c = 0 corresponding to infinitely large Reynolds number, 
can be taken giving 

V* = 2 at ; = 0 


OKlGitv.'-.v. I- : . 

05: POOR QUAUrv 
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3. RESULTS AND DISCUSSION 


The analysis has provided approximate closed-form solutions that describe the structure of 
the turbulent vortex pair generated by a lifting surface. Two regions of the wake are identified; 
a region of persistence followed by a region of decay. The initial persistence of the wake is 
associated with an acceleration of the Cow along the axis of the vortex causing a radial inflow 
which counteracts the outward radial diffusion. When the axial velocity recovers to the free 
stream value the inflow ceases and the vortex starts to decay. 

fn the interest of brevity the results discussed here will be primarily for the case of an 
el'iplically loaded w ing with brief mention of the influence of inboard and outboard loading which 
may be determined when the quantities s and e = are known. The general expressions 

for the quantities of interest are given in equations (17)-(20) which reduce to the simpler forms 
discussed below with s = q and e = J corresponding to elliptic loading. 


i lie lengiu oi persistence < 


d= io.A 
Cl 


arid shows that the length of persistence varies directly as the product of the span and the aspect 
ratio and inversely as the lift coefficient. For a typical large transport aircraft in a high lift 
configuration (fc=200ft. AH =7 and Ci=l) the length of persistence is approximately M.OOOft 
(2.05 miles). On the ether hand for a typical fighter aircraft (6=50ft, AR=l, Cg=2) the length 
of persistence is only 250ft. 


Considering now a description of the vortex, in the persistence region , the radius of the 
turbulent core sq. where the flow departs from the potential vortex is given by 


n — .1756 

independent of aspect ratio and lift coefficient. The laminar subcore radius, taken as the location 
of maximum velocity depends on the Reynolds number and (for elliptical loading) is: 



where c = *717-, Re = For Re = 10 7 ,e — .001 the value of r* is 


r* = .0186 

thus, for a transport aircraft of 200ft span, the maximum velocity occurs near the center of the 
core at a radius of 3.6ft. 


The maximum velocity is given (again for elliptic loading) by 




Thus for the transport aircraft (A/? = 7, Cg = 1) and f. r =300ft/sec the maximum velocity near 
the center of the vortex i ; 50ft/sec; for a fighter aircraft (A/? — l, C l = 2 and V ~ 300fl/sec) 
however, the maximum velocity would be 760fl/sec (ignoring compressibility effects). 


In th e decay region , i.e. for x > 10.4^6, the radius of the turbulent core and the subcore 
grows as (J)‘,i e., 

r l (r) = .173(^) i 6 = .0S4(^) { (i) } 6 


and 




For the transport aircraft given in the previous example the turbulent core extends to a 
radius of 70!t at a distance of 10.6 miles behind the aircraft and the radius r*, at which the 
maximum velocity occurs, is approximately 7ft. The velocity i j at this point is 25ft/sec. For the 
fighter aircraft the velocity rate of decay with distance is greater (since the velocity reduces by a 
factor of 2 in a distance z = ad, and d is only 250fl in this case). 
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From the general equation for d (equation 17a) and tv (equation 10a) it ean be inferred that 
configurations which have smaller values of i$ on the vortex wake (small Cl, large Aff, large a) 

also have larger distances of persistence, (and those configurations having higher values of t> have 
smaller distance of persistence). Moreover, from 10a it can be seen that in the decay region the 
only dependence on the span distribution parameter a is t> ~ a relatively weak dependence 
which suggests that the reduction of the high velocities in the far wake through tailoring of the 
span distribution of lift is not a very promising technique. 

Although comparison with experiments is not made in this paper it ha3 been found that 
the general trends are reasonably well represented by the analysis given here. The dependence 
of persistence on the axial flow at the center of the vortex probably deserves further analysis, 
however, including the influence of promoting decay through mass injection thereby permitting 
the axial flow to recover earlier to free stream Conditions. 

Finally, application of the solutions given here to rotor wakes also deserves mention. In this 
regard the effects of Reynolds number on the behavior of the velocity v $ near the axis of the 
vortex may be much more important, particularly in a proper description of the interaction of 
the vortex with a following blade. In this situation the effects of viscosity in the laminar subcore 
play an essential role in eliminating the singularity in the velocity gradient thereby permitting a 
realistic calculation of the vortex-blade interaction including the acoustic field so generated. 
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Vortex scale and velocity unchanged 
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Decay 
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Fig. I Trailing vortex flow i 
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Fig.2 Vortex structure: (a) Potential flow, (b) Turbulent flow, 
(r, = radius of turbulent core, r* = radius of laminar suheore) 






Fig.5 Radial profiles of circumferential velocity 
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INTERACTION OP A TURBULENT VORTEX 
WITH A LIP TEN G SURFACE 

D. J. Lee* and L. Robert** 

Joint Initiate for Aeronautics and Acoustics 
Department of Aeronautic* and Aeronautic* 
Stanford Univer*ity, Stanford, California 


ABSTRACT 


The impulsive noise due to blade-vortex-interaction 
is analysed in the time domain for the extreme case 
when the blade cut* through the center of the vortex 
core with the assumption* of no distortion of the vor- 
tex path or of the vortex core. An analytical turbulent 
vortex core model, described in terms of the tip aero- 
dynamic parameters, i< used and it's effect* on the un- 
steady loading and maximum acoustic pressure during 
the interaction art determined. 
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NOMENCLATURE 

speed of sound in the medium at rest 
equivalent tip area 
equivalent tip aspect ratio 
blade temichord 
nosdimeasionil constant jV r ' 
induced drag coefficient of the tip 
lift coefficient of the tip 
Theordorsea function 
characteristic length of vortex persis- 
tence 

width of the signal in time_ 
aerodynamic efficiency -- jfc - 
blade surface function or frequency Ht 
nondimtnsional circulation function 

r/r, 

Bessel function of the fint kind of order 

0, 1 

constant 0.06 or reduced frequency 
ub/U 

span length of the acoustic source 
?<fach number 
unit vector in j direction 
acoustic pressure at the position of the 
observer or pressure on the surface 
compressive stress tensor 
distance between moving source and ob- 
server or rotor radius 
H component in yi.yi.yi direction 
cylindrial coordinate 
turbulent vortex core radius 
initial r, 

laminar sub core radius 

half the distance between the equivalent 

vortex pair 

Reynolds number I*,/* 

rotor radius at maximum circulation 
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initial distance from the observer and 
source in y, direction 
nondiraenrional time Ur/b 
tip loading parameter s, = l',Giyl[R - 
*-) 

blade surface 

reception time (or observed time) 
peiiod of the acoustic signal 
correcting blade speed 
vortex velocity in r, S.z directions 
vortex velocity function 
Maximum velocity function at z = *• 
normal component of a gust velocity on 
the blade 

blade fixed coordinate in aerodynamics 
space fixed coordinate in acoustics 
initial distance from the vortex core to 
the leading edge of the following blade 
blade fixed coordinate in acoustics 
•*/r, 
r’/r t 

circulation on the blade lip or circula- 
tion in vortex core 

maximum circulation on the blade tip 
maximum circulation on the vertex = 

r_ 

delta function 

delayed reception lime due to the span 
lengtn 

eddy viscosity i’fr/p )lr 
gradient (£.£.-) 
absolute viscosity 
kinematic viscosity m/p 
air density 

air density in the medium at rest 
signal widtb/peried 

time at the leading edge during the in- 
teraction 

disturbed potential function 
angular frequency 
helicopter rotor angular speed 
Kutsner function 

emission (or retarded time) for moving 
source or circumferential shear stress 


1. INTRODUCTION 


One of the primary sources cf noise for a helicopter, 
under certain Eight conditions, arises from the interac- 
tion of the vortex produced by a blade with the aero- 
dynamics of the following blade. To analyse this impul- 
sive noise, sonrlimes called "blade slap," three parts of 







tb« phenomena must be considered; tbe tip vortex gen- 
edited by the blade, the unsteady pressure produced 
oa the following blade during the interaction, and the 
acoustic radiation due to this uniteady preuure field. 

Tbi» phenomena was first studied by Leverton and 
Taylor 1 experimentally using a rotating blade and two 
opposed airjett to simulate the tangential velocity pro- 
file of the vortex. Using the theory 2 of nonuniferra flow 
past a thin airfoil. Sears 2 treated the line vortex aj a 
gust; the vortex is forced to move parallel to the blade 
span with a displacement of half of tbe chord height. 
Parthajarathy and Karamcheti* solved this problem 
with a free point vortex, whose path is disturbed by 
the blade during the interaction from the intial height 
of half chord, and they calculated the quadrupole cl- 
iect corresponding to the acoustic sources in the Sew 
fields. Widnall,* who formulated an acoustic model for 
an oblique forced vortex using tne quasi 2-D unsteady 
aerodynamics obtained by Johnson* and Filatos, 2 ex- 
tended tbe theory to study the blade tip loading shape 
effect* and the high speed effect*. The Bets inviscid vor- 
tex model 19 was used to relate the tip loading shape and 
the viscous core was treated by using an effective dis- 
tance between the center of the vortex and the blade. A 
number of researchers 11 12 11 have recently studied tran- 
sonic effects on the unsteady aerodynamics due to the 
induced velocity cl the vortex and other! 1 * have tried 
to predict the noise cl the rotating blade using the mea- 
sured aerodynamic pressure 15 on the blade surface. 

In this paper, we attempt to study the effect cf the 
turbulent viscous vortex core on the unsteady loading 
and acoustic pressure for the case of a close encounter 
with a following blade in certain flight condition as 
shown in Fig. la. To analyse this phenomena, a simple 
analytical model of the vortex core, which is turbulent 
and viscous, is obtained by relating it to the tip aerody- 
namic parameters and we assume that the in.eraciicnj 
occur in the persistent region cf the vortex. To simplify 
the problem, we consider tbe situation where the vor- 
tex filament generated from tbe rotating tip is aligned 




Fig. 1 . Trajectories cf rotating tip vortex «’ 


parallel to the following blade a* shown in Fig. lb; in 
that situation the maximum acoustic pressure due to 
interaction it observed 1,lT . Even with the assump- 
tion of 2-D aerodynamics, the flow is too complicated 
to analyte completely as indicated in the experiment 
done by Ziada and Rockwell 1 * so we assume that dur- 
ing the interaction the vertex path and vortex core are 
not distorted when cut by the following blade. Acoustic 
pressure is calculated for the finite unsteady source mov- 
ing straight toward the fixed observer and analysed in 
the space fixed coordinate and reception time domain, 
whereas the source is described in the body fixed coor- 
dinate and emission time domain. Based on the analysis 
mentioned above, the trade-off of noise and performace 
is discussed. 


2. ANALYSIS 
2.1. ROTOR TIP VORTEX 

It is observed that the retor tip vortex quickly rolls 
into a concentrated vortex and persists for many span 
lengths before it decays (as in the case of the fixed wing 
tip vortex). However, because of the unsymmetric cir- 
culation shape toward the tip and the unrolled inboard 
vortex sheet of a rotor blade, the rotor vortex is more 
complicated to describe. Thus, for simplicity the ro- 
tor vortex field is modelled as a vortex pair, i.e., a tip 
vortex with an equivalent counter vortex of the tame 
strength , which replaces the unroUod inboard vortex 
sheet. The earratute of the vortex filament is neglected 
locally during the formation and the oncoming veloc- 
ity near the tip is assumed equal to the velocity at 
the position cf tbe maximum circulation oa the span. 
With this tip vortex model, the vortex core structure 
and strength are related with the tip aerodynamic pa- 
rameters and are determined by following the general 
approach used by Spreiter- Sacks 1 * for the potential part 
and by Roberts, 19 , for the turbulent viscous core. 

Relation between tip aerodynamic parameters 
and vortex core 

Consider a rotor blade cf radius, A, rotating with an 
angular velocity R in a stream cf uniform velocity U 
The blade will have a point of maximum circulation 
at a radial distance A*. With the assumption cf an 
equivalent tip vortex, it is convenient to consider the 
lift and drag on the rotor outboard of the radius A.. 
Thus 

(1) 

Generally, the lift coefficient Cl, and equivalent lip area 
A, are the functions of asimuthal angle v. Tbe lift can 
be expressed, alternatively, as the integral of the distri- 
bution of circulation around the blade ax 

L,„ = p[f»R,) • 2(A - Rm)C m £ C) 

where F_ is tbe circulation at A,, and C * r/r„ is 
shown in Fig. 2. Then, from (1) and (2), 

T- w(nA„){A-A.)^-i (3) 
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Fig. 2. Relation between the spanwise circulation 
and sip vortex 


where 


7> ( Fr 


Rm 


w 


Tbt quantity », represents the load characteristic cf the 
tip and A R is the equivalent aspect ratio with equiva- 
lent span length of 2 (R - £_)• 

It is well known that a lise cf the rolled up vortex is 
related to the lift-induced drag and is determined by the 
integration, in the plane normal to the axis c,f the vor- 
tex, of the rate of formation cf rotational kinetic energy, 
as wu nted by Sprtiter 1 * t.e. 


D " m !L.W + vi’)iS + 2 j ' lpv]2rrdr 


where 

and 


d„ - gp(n&.) , c„, 1 / 1 
* 


V# ** 


_r_ 

2xr 


( 5 ) 

( 6 ) 
(7) 


The rotational energy has two parts; one from the po- 
tential part of the two point vorticies outside cf the 
equivalent tip vortices and the other from the vertex 
cores of radius r,. The vertex cere strucurs is assumed 
not affected by the vortex pair, so it becomes a circu- 
lar £ce and its circulation T is a function of the ra- 
dial distance r. The strength of the point vortex, T| is 
assumed the same as the maximum circulation on the 
blade, r.,which is actually almost CG95 of r«, in the 
experiment data of Tung”. 

By combining Equation (2) and (7), Equation (5) 
can be expressed as below after integration cf the po- 
tential part. 


? '£) t * 

where the quantity c, is the tip aerodynamic efficient 
defined is 

Ct? 


' ».<AC 0 „ 1 

and r t is half of the reparation distance between tl 
vortex pair determined by using the conservation of ce 
troid cf the vcrtici.y 


'l D (S- ft-) f t = (/J - R.)s, 


(10) 


Then r, f the rvdiu* of the turbulent vortex core, nay 
be expreued in cloied form u ikewa below 


r, w (J? - £«!». |sinh - j|jj (11) 

The rotational energy cf the core is obtained from the 
G in equation (16) , which will be discussed in the fol- 
lowing section. 


The turbulent viscous vortex cor* 

For large distances from the vortex center, the tip 
vortf.': behaves like a point vortex cf constant strength 
T|. fltvr the center of the vortex, where the turbulent 
and viscous effects are significant, the vortex cere can 
be described by a differential equation expressing the 
transport of angular momentum, 2erv#, in body fixed 
cylindrical coordinate as shown in Fig.3. As- 

suming the axial gradient is smaller than the radial gra- 
dient, the equation is written 


SG SG 9 It . t 13C ! 

u — — 4- v— = r— l(v -a e)-— (12; 

dx dr dr 1 ' 1 r dr 


where u is the axial velocity and v is the radial velocity, 
and t \c are the kinematic and eddy viscosity, respec- 
tively. The eddy viscosity for vortical flow has been 
aaalysed by Hoffman and Joubert” in a way aaalagous 
to the traditional muting length theory of Prandtl : 

c~k\ r/ P )\r (13) 

where 

*/e “ (><) 

(k » 0-06 appears to fit the experimactal data” at a 
distance, d where the vortex starts decaying). 

In the axial direction, there it an axial velocity de- 
fect in the core due to the viscosity such that the mean 
velocity in the plane normal to the axial direction is 
zero at the trailing edge. It is assumed here that the 
axial velocity viriei linearly from the trailing edge to a 
distance, d, where the vortex starts decaying. Then, the 
inward radial velocity v is determined from the conti- 
nuity equation. This inflow into the center balances the 
viscous and turbulent diffusion such that a persistent 
region of the tip vortex exists. 

With the velocity components v,v assumed and 
eddy viscosity formulated in equation (13), a similar- 
ity solution far G was obtained by Roberts 1 ’ with the 
boundary conditions that the tangential velocity, »> at 
the center is zero and the core merges with the poten- 
tial flow. The tangential velocity of the vorttx can be 
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Fig. 4. Vortex velocity distribution 
expretsed u follow! : 


— isfei 1 ™ " s > 

where 

V(i) o C(*)/* 

- [l -(I -I*) (l + 7)']/* 4<t (1C) 

■ 1/a 1 > t 

The similarity variabls, a. i* defined u r/r,(n), where 
the turbulent core radius r, remains constant and equal 
to r,(0) for r < d and increases as r,(0)\/r/d for r > 
d. The persistence distance, d is x/«* : (0/l.r;(0)/r,). 
The velocity function V(i) of equation (16) ii valid for 
small c .which ij a ncndimcasionil parameter related to 
the inverse to the Reynold! Number, Re, denned aa 


7ru 2r J_ 
*T, “ Rt 


(17) 


where Rt, baled on the maximum circulation in the 
vortex, is related to the blade chord length as a charac- 
teristic length : 


2.2. UNSTEADY AERODYNAMICS 
DURING THE INTERACTION 

’Vith the tip vortex defined ia the previous section, 
the unsteady flow field during the interaction between 
the vortex core and the following blade will be described 
with the assumption of no distortion of the vortex path 
or of the vortex core. Even though the vortex is turbu- 
lent and viscous, it is assumed that the flow around the 
airfoil is a potential flow during the interaction. These 
assumptions mean that the vortex is forced, not tree, 
and the vcrticity field of the vortex core and the int>- 
tational field around the blade can be split and cou- 
pled through the surface boundary of the blade aa ex- 
plained more generally by Kovissnay 1 ’ and Coldstein 1 *; 
this ia the tame aasumption as for traditional gust the- 
ory. Thus tbs vortex behaves likt a convecting wave. 
These concepts were visualised by Rcgler 5 * ia the prob- 
lem of the interaction of a vortex array, representing 
the turbulent flow, with a semi-infinite plate. It shows 
that tht stream function representing the vortical fields 
is not distorted but the combined stream function in- 
cluding the disturbed irrotalional fields, satisfying the 
boundary condition of no Sow through the solid bound- 
ary, ia distorted to give a vortex on either s.de cf the 
plate. This flow pattern is in genaral agreement with 
the experiments cf Yu 1 * in the extreme case for which 
tha center of the vortex cor* meets the blade. Within 
thest asaumpliens, only the normal component cf the 
vortex velocity cn the blade surface is important espe- 
cially ia case of no angle cf attack between the conven- 
ing flow and blade. Even though the assumptions art 
valid only when the normal component cf the velocity 
is small compared to the convening velocity, we use the 
assumptions generally to dcrivs the unsteady loading on 
the blade surface during the interaction. 

Unateady loading on tha ear-face 

The governing equation for the incompressible small 
perturbation potential fuclion, 4, is 


R * “ — (i6) 

k AR, 1 , v * ' 

(for a typical model rotor blade e -< 0.01 to 0.1 and 
Rt - 10* to 10‘). 

The velocity fuction Vfx) can be expressed in a sim- 
ple form as shown in Fig. 4 . : 

V(*) m -j for luniajf , u b core 
y(i)»l-j for the turbulent core 

The maximum velocity V'- » 2 — 2x* occurs at the po- 
sition a" m (c/2 logl/e) 1 ' 5 as shown ia Fig. 4 and the 
velocity profile near the center, in the laminar sub core, 
shows that the velocity gradient is steeper for smaller 
value cf c; that is, for large Rt, which property becomrs 
important for close encounters of the vortex with a fol- 
lowing blade. 


V*d w 0 

and the boundary conditions are 

~ w -ic,(l - Ur) |r|<6 l»0 

OX 

p w 0 | r |> 6 1 «= 0 


(20) 


(2») 



Fig. 5. Coordinate system of the blade during interac- 
tion 
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where u, i« the normal component of the vortex fust. 
The preisure is given by the Bernoulli equities 


(34 „S6\ 

* u r IJ 


( 22 ) 


Using the conservation of circulation in the Cow field 
ltd the Kutta condition it the trading edge, the lift for 
a tharp-edged gust w u obtained by Kussner \ (with the 
assumption that the shedding votticies from the trailing 
edge covect wi:h the nine velocity of the free stream), 

I(s) - (23) 

where m constant and 


<(») 

and 


2v J-m •« 


j.(*re 


•tit. is 


-dk 

(24) 


s 


Ur 
b ' 



(25) 


The variable s is ccndimeasiosaliied time and k it the 
reduced frequency. The (unction c ( * ) is the bf: response 
function for the step gust, called the Kustner function. 


For the vortex gust, we can use the Duhamel theo- 
rem to obtain the lift using the Kiissner function. Thus, 


tr.(0)ti(s) -e)io 


I(s) » 2*e.Vb 

[ 49 tW 

(20) 

Here, the velocity of the vortex gust u >, is obtained after 
the transformation of the coordinates as shows in Fig.S 


and the gradient of the gust velocity is 


(27) 


dw,;?; ^ t v , i. - ot> 

dj 2rr , f| r, 


(2S) 


where the quantity at time <t is of the converting vortex 
velocity at the leading edge and r. is the initial distance 
between the center of the vertex and the leading edge. 
As ib own in Fig.6. the gradient of the gjst velocity de- 
pends on the turbulent vortex site, r, t the lam' \r sub 
cere radius, r* acd the slope at the center of the c v 2/c. 
The lift can be evaluated easily analytically or numeri- 
cally if the Kustner function, which has the integration 
for the whole range of frequency as shown in equation 
(24), can be appreximcted in a simple analytical form. 


Behavior of the KiUsner function 

The approximate form of the Kiissner function can 
be expressed in a series of exponential terms as '* it 

v(s) « l-0.Se- ,,u -0.Se- s>2 (23) 


This is valid for the whole range of s except near the 
origin where a mare accurate form can be obtained by 
using the aaymtotic expression r! the function 1 near the 
origin, s < 2, (the time required for the front of the step 
gust to pass ever the chord), as follows 


<(*) 



£ 
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s < 2 (30) 


Comparing (29) and (30), there is no significant differ- 
ence in the unsteady loading as shown m Fig. 7 . But 
if this loading acts as the acoustic source, a vgmficant 
difference occurs near the origin because it involves the 
gradient of the function. Then for s < 2, 


difs) t\\ 3 J_ 

ds “ r Us” 2 12 


S I 
2 90* 


7 23 
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and for s > 2, 


- 0.5 0.13«-* ,5 '"‘ 1 " ~0.Se-‘~' 5 " (32) 

as 

The gradient of the function has an integrable singular- 
ity at the origin and the exponential terns art modified 
to give a smooth function at s = 2. This is shown in 
Fig. 8 . 



Fig. 7. Kustner function 



Fig. 6. Profile of the coavecting vertex velocity gradient 
at the leading edge (cemtr thif’.ed by x,/6) 


d*(s) 

dt 
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Fig. 8. Gradient of the Kussner function 








Fig. 9. Coordinate of th* moving source 


2.3. ACOUSTIC PRESSURE 

The acoustic pressure field for a moving body in 
the pretence of a vortex can be obtained by tiling 
the Ffowcs-WiUiams and Iliwkings equation 51 in the 
medium at reit and by assuming the nonunifcrmity ii 
confined to the sourcei. The aignal is aealyted in space 
fixed coordinates and the sources are described in body 
fixed coordinates. The equation has monopole, dipole 
and quadrupol* sources corresponding to thickness, sur- 
face pressure and unsteady Jos fields. Here we consider 
only the dipole source because, assuming a thin airfoil 
and low math number, we can neglect the moaopole and 
quarupole terms. As explained by Curie”, the dipole 
source in the equation implies the effect of scattering 
due to the body for wave lengths of the signal large 
compared with the body site; then the source can be 
represented by the total force on the surface. Actually 
the observed acoustic signature during the blade vor- 
tex interaction contains all frequencies’ 1 but considering 
that most of the energy comes from the lower frequen- 
cies, a compact source can give the basic property of 
the signal. Here, in our problem of a finite source with 
large AR, a spanwise noncompact source is assumed. 

The inbomogenous wave equation with a moving 
dipole source is 



where a. is the speed cf sound in the medium at rest 
and p,, is the compressed stress, and / is a function of 
the body surface. 

The solution in 3-D space is 


r(x.O 


1 3 r P.,», 

irdz, h [/} | i _ m . i | 


dS(y) (34) 


where p = aj(p-p.) and R =| x-tl/r, - y | as shown in 
Fig 9. The source is evaluated in emission time r,. The 
time relation 


is related to the Doppler shift because of the moving 



source and 


(36) 


1/ I 1 - Af ■ R/R | (36) 

is the Doppler factor in amplitude. 

Far field approximation 

For the far field , the solution can be expressed in 
terms of a time derivative instead of a space derivative 
and for a chordwise compact but spanwise noncompact 
source, it can be written as follows (after using the re- 
lation cf (35)) 

! ft? 


p(*.0 = - 


<*a. fi-'ll-M-J ?; 1 

C [£*('. vs)] /y. 


where R' is the distance from the observer to the cen- 
ter of the source. For rectilinear moving source, the 
emission time it obtained as shown below 

t - 

j — .. — *t. 

* 1 - Af» 


/(« - ^) ! - o - m {«* 

1 1 - so 




The acoustic pressure is composed of two parts; the di- 
rectivity function of the point source and integration of 
the source along the span. The integrand of the lime 
variation of the lift can be obtained at shewn below 
using the Duharael integration of the gradient cf the 
Kussner function, defined in (31) and (32), which rep- 
resents the acoustic response function for the step gust, 
dl 

^(s.y,)=2irp.Ui 




yi ) di(s - 

di 


The acoustic pressure also depends on the gradient of 
the gust which is function of r, and c as shown in Fig.6. 
The schematic diagram for the Duhamel integration is 
shown in Fig.10, where the contribution of the gradi- 
ent of the center of the core is large because of the 
de'ta-function-like acoustic response fucticn. The con- 
tribution to the integral along the span of equation (37) 
comes from two factors; one from the lift variation along 
the span, which is sero in our 2-D aerodynamics, and 
the other from the emission time variation along the 
span as seen in equation (33). 
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Fig. 10. Schematic diagram cf the Duhamul’s 
integration 
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An alternative form of the integral in equation (37) 
cm be obtained by changing the variable of integration 
Vi to r u «hown below 


**- t} ~ / Tt L{,) 


i'dLlljr 

3r 


(«) 


Then, from the relation of the variable! y s and r in (38), 


__ dLU) {(1 

dt Afv/Jt-J)’ - Ar((«;-») : -^ ! ) 

Ml) 

where the quantities 7JT and 757 are tbe noadimension- 
alised distances by the semichord length b and 1 is the 
noditnensionalised reception time by b/U at s in (25). 
This form is useful in evaluation of the integral because 
the integTand is the direct function of the variable of 
integration, s. 


Interference of the ipanwice non-compact lonrce 

The integral of the spanwise noncompact source 
can be interpreted as the interference between multiple 
point source signals. Then dL(t)/de can be expressed 
in the emission time domain with measured or com- 
puted values of Z(s) for each span location, Each signal 
is transformed to the reception time domain as shown 
in Fig.l 1(6) with the time relation t = r, £(r,,y s )/a.. 
The amplitude is the same for corresponding f and r. 
Then each signal is summed in the reception time do- 
main as shewn in Fig. 11(c) . The delayed lime (in the 
reception time domain) due to the different source po- 
sitions causes the interference, which contributes to the 
amplitude in a destructive way and also to the directiv- 
ity pattern in adddion to the point source directivity. 

The interference is closely related to the ratio of At 
and d, (nondimensionalited by b/U), where the width 
of the signal, d, depends on the position of the maxi- 
mum velocity of the vortex and the maximum delayed 
time along the span, At depends on the observer posi- 
tion, as shown in Fig. 12. If the observer it far from the 
plane of symmetry of the source, greater interference is 




Fig. 11. Interference process of the spanwise noncompact 
source 



Fig. 12. Variation of delayed signal time with observer 
position 




Fig. 13. Simplified model signal in time domain and its 
spectrum 


maximum acoustic pressure is observed; then the source 
can be approximated as a point source to evaluate the 
maximum acoustic pressure in this simple geometry . 

The effective ratio of Af and d, can be estimated 
from a simplified point source signal and its spec- 
trum as shown in Fig. 13. Assuming that the spec- 
trum begins to decay significantly at np, a N where 
the the solidity p, is defined as d,/T ,, the correspond- 
ing angular frequency of nth(= ,%'//>,) harmonic, u„ is 
2 t n/(T,b/U) = 2rN/[d,b/U). Thus the reduced fre- 
quency, *„ is u.b/U =c 2 r.V/d, and the frequency, /„ 
is u/„/2x = ,V/(d,6/b’), and the nondimensionalited pe- 
riod, T, is l//„ U jb :<■ d,/.V. From Fig. 13, it is seen 
that for N = 1, then T, «= d,, the contributions of the 
waves for which periods are less than d, are small. So if 
the delayed time is greater than the period T« (approx- 
imately, Af > d,/4), interference will occur. 


3. RESULTS ic DISCUSSIONS 

From the previous analysis, we find that the un- 
steady lift and the acoustic pressure are closely related 
to the slope of the vortex velocity, which is the function 
of the vortex turbulent core si:e r, and a nondime nlional 
parameter c (inverse of Reynolds Number Re). 

The vortex velocity profile is computed and com- 
pared with tbe experimental data obtained by Tung”. 
The unsteady lift and acoustic pressure are calcu- 
lated and the variation of maximum acoustic pres- 
sure with the parameters r, and c is examined. Tbe 
coise/performance trade-off for blade tip loading shape 
is also discussed. 
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Vortex Velocity Profile 

Tte vortex structure generated by a blade having 
radius R = 1,05m and aspect ratio AR = 13.7, are 
computed. The general structure shows good agree- 
ment with experiment even though the maximum ve- 
locity and it> position are underpredicted aa shown in 
Fig. 14. The velocity profile predicted with our analy- 
sis is compared with that predicted w,th Bet: analysis 
lor the persistence region. Both show generally good 
agreement outside of the core but within the core the 
present analysis is more realistic since the Bets vortex 
has infinite velocity especially for the elliptic tip load- 
ing. In both cases, the distance from the location of 
the maximum circulation to the tip is talcen to be 5% 
of the blade span length. It should be mentioned that 
the maximum rotational velocity of the vortex is almost 
40% of the rotating tip speed in this cane. But if the 
position of maximum circulation is shifted toward the 
root, the maximum velocity of the vortex is reduced 
substantially. 

Unstesdy Lift and Acoustic Signal 


w,(s. - s») 


do. 


(*. - st) 



Typical nondimensional unsteady lift and acoustic 
pressure for the vortex gust are plotted in Fig. 15. It 


Fig. 15. Typical variation of unsteady lift and acoustic 
signal 
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(a) velocity profile 



Fig. 14. Vortex core velocity profile comparinon 


is observed that the lift profile is similar to the vor- 
tex velocity profile but distorted due to the behavior of 
the Kuasner function which introduces a time lag near 
the origin. The acoustic signal, in the tmiasion time 
domain, is similar to the slope of the vertex velocity 
profile because the acoustic function acts like a delta 
function which assumes the value of the slope during 
the integration. So the amplitude of the peak acoustic 
pressure is related to the slope of the vortex at the cen- 
ter of the laminar sub core, 2/c and the width of the 
signal is relate- to the radius of the laminar core, r*. 

Reynolds Number and Vortex SIse Effect 


It is to be expected that both the peak unsteady lift 
and the acoustic pressure decrease with larger turbulent 
vortex core site r, and smaller Reynolds .\umber(/fe = 
Ti/v) because the gust function is r,/Irr, ■ Vfir/r, ). 
But the acoustic pressure in more sensitive to the slope 
of the core(/Ie — 2/c) showj in Fig.6, because again the 
behavior of the delta function like the acoustic function 
varies as the slope of the laminar core, 2/c. For the 
maximum acoustic pressure, the effects of tbs Reynolds 
Number and vortex siie are shown in Fig. 16. To re- 
duce this peak pressure, it it obvious that we need to 
reduce the vortex ei:e r, and Reynolds number Re. One 
method would be to increase tbe distance, [R - 71.) be- 
tween the location of tbe maximum circulation and tip, 
because r, and Re both are directly related to it from 
the equation (II) and (13). But this would result in a 
Ions of lift near the tip which is proportial to (n/?„ )* 
with tbe same lift coefficient. The other method would 
be to change the tip loading shape which is related to 
r, for the fixed value of {!} - /L.) and T t . 


[ 
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Fig. 16. Effects of Reynolds Number fit and vortex 
site r, t os maximum acoustic pressure 


Noil* Performance Trade-Off 

The influence of blade loading on the noiit and per- 
formance bai been itudicd. The variations of the cir- 
culation ibape of the form r/r, = (sinS) 1 "" 1 , where 
cos 8 = y/(R - R«), are considered as ihown in 
Fig.l8(a)(m = l is elliptic loading). For large value 
of m, the core radiui increases and there is a reduc- 
tion of noise compared with the elliptic loading, but 
the tip loading relief also gives a reduction in the tip 
aerodynamic efficiency defined in Equation (9) as shown 
in Fig. 18(6). From Fig. 19, a maximum noise reduction 
S iB is obtained for m = 2 with 25% loss of efficiency 
and a noise reduction of 3d£? with 5% less of efficiency. 
For the triangular loading there is more than 8 iB re- 
duction in noise but also more than 30% reduction in 
efficiency. 



(e) 



Fig. 17. (a) Tip circulation profile (b) variation of f\!R- 
&m and &p(dB) with tip-loading parameter ra 
(c) Variation of tip efficiency with tip-loading 
parameter m 



Fig. 18. Trad* off between nois* reduction and aerody- 
namic efficiency for various m. 


4. CONCLUSIONS 

The blade vortex interaction where the turbulent 
viscous core is cut through by the blade is analysed with 
the assumptions of no distortion of the vortex path and 
or vorticity of the core during the interaction. From the 
analysis and results, the following conclusions can be 
drawn : 

(1) The detailed vortex structure including the viscous 
and turbulent core generated by the rotor tip can be 
approximately predicted with a simple turbulent core 
model and equivalent tip model. 

(2) The interaction between the vortex core and the 
blade is analysed with the gust analogy under the as- 
sumptions of no distortions of the vortex for the extreme 
case where the center of the core is cut through the vor- 
tex. 

(3) Acoustic pressure is obtained by introducing tbe 
acoustic response function for the step gust for the 
chordwise compact source. It shows that the ipanwise 
non-compact source gives a sisable effect when the ob- 
server is off-center of tbe source because of the Urge 
interference effect. 

(4) Maximum acoustic pressure varies with Reynolds 
Number and vortex sire, which are related to the tip 
loading. 

(5) Shape modification gives a maximum reduction of 
6 ~ 8 iB including the triangular loading, but a reduc- 
tion of 3 JI3 in maximum acoustic pressure compared 
with the elliptic loading, is expected for 5% loss in effi- 
ciency. 

For the extreme care of cutting through the core by 
the blade, the chordwise non-compact source should be 
considered in order to obtain results over tbe range of 
frequencies in the high speed case. It should mentioned 
that the non-linearity due to the large maximum vor- 
tex velocity and due to the quadruple for this case, 
which would be the same order of tbe dipole, should be 
considered also. 
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ABSTRACT 

An attempt has been made to model the so called “Leading Edge Vortex" which 
forms over the leading edge of delta wings at high angles of attack. 

A simplified model has been considered, namely that of a two-dimensional, 
inviscid, incompressible steady flow around a flat plate at an angle of attack with a 
stationary vortex detached on top, as well as a sink to simulate the strong spanwise 
flow. 

The results appear to agree qualitatively with experiments. A comparison has 
j also been made between the lift and the drag of this model and the corresponding 

j results for two classical solutions: 

j 

(i) that of totally attached flow over the plate with the Kutta condition satisfied at 
| the trailing edge only. 

j (ii) the Helmholtz solution of totally separated flow over the plate. 

j 

| 




J 
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• 
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NOMENCLATURE 


/ 

/ / * 

English letter symbols: 


• 

a 

radius of the cylinder. 



b 

length of the span of the plate. 



c 

chord length of the plate. 



Co 

drag coefficient. 



C L 

lift coefficient. 



Cm 

pitching moment coefficient. 



c P 

pressure coefficient. 



d 

distance between the streamlines leading to the two stagnation points on 
the plate. 



D 

drag force. 

; 


F 

total force on the plate. 

i 

i 


k 

total vortex strength. 

! 


*0 

, bound vortex strength. 

! 



leading edge vortex strength. 

» 

j 


K 

non-dimensional total vortex strength. 

J 


Ko 

non-dimensional bound vortex strength. 

j 


Ki 

non-dimensional leading edge vortex strength. 

i 


L 

lift force. 

| 


m 

mass Qow rate. 



m, 

sink strength. 



M, 

non-dimensional sink strength. 



n 

non-dimensional parameter showing position along the plate. 



r,0 

porlar coordinates in the z-planc. 



Ri 

non-dimensional distance from the center of the cylinder to equilibrium 
point. 



5 

wing area. 



U, V 

velocity components. 



Voo 

free stream velocity. 


l 

w 

complex potential. 



x,y 

cartesian coordinates in the z-planc. 



k> 

complex variable in the original circle plane. 


/ 

** 

complex conjugate of z. 



z' 

complex variable in the rotated circle plane. 


/ 

/ 
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Greek letter symbols: 
a angle of attack. 

T total circulation. 

To circulation due to the bound vortex, 

f complex variable in the final plate plane, 

f' complex variable in the Joukowski plane. 

coordinates in the ^-plane. 
f'.q' coordinates in the f'-plane. 

Ao doublet strength. 

p air density. 

Subscripts: 

e in the plane of the cylinder, 

p in the plane of the plate, 

r radial component, 

sc on the surface of the cylinder. 

»p on the surface of the plate. 

TE at the trailing edge. 

z in the z- plane. 

0 tangential component. 

f in the f-plane. j 

0 refers to the origin {center of the cylinder). 

1 refers to the equilibrium point. 

( 

♦ 




TABLE OF CONTENTS 


Abstract i 

Nomenclature . ii 

Table of Contents iv 

List of Figures vi 

1. Introduction 1 

2. Analysis of Flowficld . . 2 

2.1 The Complex Potential 2 

2.2 The Velocity Field 3 

2.2.1 Radial Components 3 

2.2.2 Tangential Components 3 

2.3 Conformal Transformation Used to Analyze Flowficld 3 

3. Conditions 5 

3.1 Kutta Condition at Trailing Edge 5 

3.2 Kutta Condition at Leading Edge 5 

3.3 Vortex/Sink Velocity Condition 6 

3.4 Non-Dimensional Parameters 8 

3.5 Additional Condition 9 


& 


53 



Results 11 

4.1 Position of Equilibrium Point 11 

4.2 Vortex and Sink Strengths 11 

4.3 Pressure Distribution on the Plate 11 

4.4 Lift, Drag, and Pitching Moment 13 

Conclusions 15 

References 15 

Figures 17 

Appendix 1 36 

Appendix 2 37 


v 


3U 



LIST OF FIGURES 


Figure 1. Schematic of the flow over a slender sharp-edged wing 18 

Figure 2. Steps in mapping sequence from original circle plane to final flat 

airfoil plane 19 

Figure S. Velocity components in cartesian and polar coordinate systems. . 20 

Figure 4. Force diagram on the flat plate 21 

Figure 5. Streamline pattern around the fiat plate at an angle of attack with 

a detached vortex and a sink 21 

Figure 6. Position of the equilibrium point (along the plate) 22 

Figure 7. Position of the equilibrium point (perpendicular to the plate). . . 23 

Figure 8. Leading edge and bound vortex strength versus angle of attack. . 21 

Figure 9. Sink strength versus angle of attack 25 

Figure 10. Total vortex strength versus angle of attack 26 

Figure 11. Pressure distribution along the plate 27 

Figure 12. Schematic of the flow over a flat plate at an angle of attack, 

according to three different models 31 

Figure 19. Lift characteristics 32 


35 



Figure 14- Drag characteristics 33 

Figure 15. Lift to drag ratio versus angle of attack 34 

Figure 16. Pitching moment charaxteristics 35 


^’■si 



1. INTRODUCTION 
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ikwmnvv.. 


v*:a-wiv , >v' 


So far there are two well-known models for the Sow over a Bat plate at an angle 
of attack. That of totally attached flow with the Kutta condition satisfied at the trailing 
edge only, and that of totally separated Cow (Helmholtz solution). 

The present model, considering partially separated Sow lies somewhere between the 
two and despite the fact that it too is two-dimensional, gives some very useful represen- 
tation of the three-dimensional Cow over delta wings. On such wings the leading edge is 
usually sharp, causing thus Cow separations even at moderate angles of attack. These Bow 
separations take the form of two free vortex layers joined to the leading edge of the wing 
and rolling up to form spiral shaped vertices above the upper surface cf the wing (Figure 
la). 

These vortices induce additional velocities at the upper surface of the wing. The 
corresponding pressure distribution shows distinctly marked minima beneath the vortex 
axes (Figure lb). Accordingly, an additional lift force occurs which depends nonlinearly 
on the angle of attack (Figure lc). 

In the theoretical study made here, a simpliGed model has been considered, namely 
that of a two-dimensional, inviscid, incompressible steady Dow over a Cat plate at an angle 
of attack, with a stationary vortex detached on top as well as a sink to simulate the strong 
spanwise Cow caused by the pressure gradient due to sweep in the three-dimensional case. 
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2. ANALYSIS OF THE FLOW FIELD 

In the original circle plane (also called z-planc) the flow field consists of the following 
components: 

(i) Uniform wind Voo coming from the negative i-axis. 

(Li) A doublet Ao at the origin to simulate a circular cylinder j^| = a. 

(iii) A bound vortex &o at the origin to account for the circulation To around the plate. 
Note that although the flow is aligned with the x-axi3 the plate is at an angle of 
attack, requiring thus circulation To for the Kutta condition to be satisGed (sec 
Figure 2a). 

(iv) A steady vortex fcj of fiuite radius placed at zi(ri,#i) to simulate the leading edge 
vortex. 

(v) A steady vortex — k t at the inverse square point of zj induced by the circle theorem 
(see Appendix Al.l). 

(vi) A steady vortex ki at the origin also induced by the circle theorem (see Appendix 
Al.l). 

(vii) A sink mi (mi < 0) placed at z ( [r lt 0i) to simulate the spanwise flow along the 
center of the vortex. 

(viii) A sink mi at. the inverse square point of zi, induced by the circle theorem (sec 
Appendix A1.2). 

(ix) A source —mi at the origin also induced by the circle theorem (sec Appendix A1.2). 
2.1 The Complex Potential 

For the components described above (regrouping them together), the complex po- 
tential is given below: 

VJ = ^z + — ^ +(— m, + \k) hr z-f (mi+i/:i) ln(z- zi) + (m t -iki) In ^z - — ^ (1) 
where 

k- k 0 + ki (2) 


2 
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2.2 The velocity Field 

Differentiation of the complex potential gives the velocity field: 

~ = (n r - t'uo) e~'° (3) 

which applied to equation (1) gives: 


2.2.1 Radial Components 




m i[ r — r i cos (0 — <?i)] + fct r i sin(0 — 0 j) 


r 2 + rj — 2 rri cos (0 — 0 j) 

m ir*i frrj - a 2 cos($ - <?i)j — fcjrja 3 sin(0 — 0i) 
r 3 r 2 a* — 2 rr A a 2 cos(0 — 0j) 


( 4 ) 


2.2.2 Tangential Components 

„ k . mirisin(0-ffi)-fci[r-r,cos(ff-(»i)] 

tiff — — . oo l t + — 1 Sin v — r ~ s r-r -r-r 

\ r 2 / r r 2 + r J — 2rrj cos(0 — Qi) 

miria 2 sin(t? - <?i) + KirJrrj - a 1 cos(0 — 0 t )] 


( 5 ) 


r 2 rj + a 4 - 2rria 2 cos(0 — 0 1 ) 


2.3 Conformal Transformation Used to Analyze the Flow Field 

A solution for the How field is provided by a mapping sequence tiiat transforms the 
potential flow about r. circle into a flow about a flat platn at an angle of attack with a 
detached vortex/sink combination. The steps are the following (sec also Figure 2): 

a. original circle plane (z-plane), 

b. rotated circle plane (z'-pianc), z' = ze ,a , 

c. Joukowski plane (f'-planc), j' = z' + 

d. final plate plane (y-plane), f = f'e -10 . 

The first step is a simple rotation in order to make the flat plate depicted inside the circle 
aligned with the i-axis. The second step is the Joukowski transformation which transforms 
the circle into a flat plate. The third step is another rotation which gives the plate an 
angle of attack with respect to the horizontal free stream velocity. 
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Combining the three steps wc have that 


3 


and if wo set 


we get 


*=§ 

s = $ + 'n 


£ = r cost? 4- — co3(0 4- 2a) 
t) = r sin 9 sin(0 + 2a) 

On the surface of the cylinder (corresponding to the surface of the plate) r = oso 

£ ap = a [cos 0 -f cos ( 0 4- 2a)J 
rj ay = a [sin 0 - sin [0 + 2a)j 


now 


and 


dw dw dz 

df dz df 

l 


( 6 ) 

( 7 ) 

( 8 ) 

( 9 ) 

( 10 ) 

( 11 ) 

( 12 ) 

(13) 


dz _ 

df df/dz 

z 7 e 7ia 

“ z 2 e 2 *“ - a 3 

r 2 e i(2 a -h2C) 

= ~3gip a + 2P) _ a 3 

Considering the surface again (r — a), expanding the right-hand side into sines and cosines 
and using equation (3) we get Snally 


T.. 1 _ ain (‘10 -f da) — 2 sin (20 4- 2a) — 2 cos (20 4- 2a) r , 

2 4- 2 cos (40 + 4a)-4 cos (20 4- 2a) 

K.J; = ~ 2 cos (20 4- 2a) 


(14) 

(15) 
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3. CONDITIONS 


Boundary conditions in the final plate plane require that: 

(i) the flow depart smoothly from the trailing edge; 

(ii) the flow depart smoothly from the leading edge, i.e., the Kutta condition must 
be satisfied at both edges of the plate. 

Also, the vortex/sink combination must be located at an equilibrium point in the 
flow field, that is, a point where the velocity induced by all other components is zero. 

3.1 Kutta Condition at the Trailing Edge 

From equations (14) and (15) with 0 = —a we find the velocity components at the 
trailing edge of the plate: 

K.r = -5K.l‘.-~ ' (.*•) 

k£Mk.i:-- (in 

For the Kutta condition to be satisfied there (i.e., u r finite and txo = 0) we must have 

K..r-°=o (is) 

or, from equation (5) with r = a and 0 — —a 
(‘iFooasin a - k) [r 2 + a 2 - 2r 1 aco3((7 1 + a)] - 2miria sin(0i + a) -f &i(rj - a 2 ) = 0 (ID) 

3.2 Kutta Condition at the Leading Edgo 
Similarly, for the leading edge wc must have 

{uo..\°r*- n =0 (20) 

and again from equation (5) with r = a and 0 = tt — a 

(2 V m a sin a + k) [r 2 + a 2 + 2riacos(0! -t- a)] — 2miria.sin(0i + a) — &i(r 2 — a 2 ) = 0 (21) 


hi 
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3.3 Vortex/Sink Velocity Condition 

Equation (3) cannot be used to find the conditions which make the velocity vanish 
at the location of the sink and vortex because the velocity given by equation (3) is singular 
there. Threfore, the usual limiting process has to be used to find the velocity components 
at the center of the vortex/sink combination. The strengths of the singularities can then 
be adjusted so that they are stationary in the presence of the plate. Following the analysis 
of reference 3 the complex velocity function in the final plate plane at the equilibrium 
point b found to be 



-{mi+iki) 


\d 3 z/d^ 3 } ( , 

[ d */ d s](, 


= 0 


( 22 ) 


The derivative ^ has been found in Section 2.3: 


dz _ (ze' a ) 3 

rff (ze*°) 2 — a 2 


while 

d 3 : _ d / dz\ _ d T z 3 e 3,a 1 
df 2 d$ \df / df [z 2 e 2,a — a 2 j 

and substituting z from equation (7), taking the derivative and substituting back f from 
equation (6) we get 

<Pz a 2 e 2 "»[l;r 4 4- 3 (z 3 + a 2 C - 2, “) 2 ] 

• d; 3 ~ 2z (a 2 e- 2iu - z 2 )(z 2 e 2 "* - a 2 ) 2 

Also, from Figure 3 wc have the following transformations: 


ti r = u e cos 0 4- v c sin 0 
u 0 = v e cos 0 - u c sin 9 
or 

. ti e = u r co3 0 — ug sin 0 

v c = u r sin 0 4- ug cos 0 

Substituting all these into equation (22) we get 


[(u,i cos — ugi sin ) — * (u r i sin 0i 4- uoi cos0i)| ZiC 2,a 

r _ . Q 2 [•»*? -4- 3(z 2 4- a 2 e~ 2,u ) 2 ] _ 

( i + ^ 2ei(r 2 — a 2 e _2,a ) 2 ~° 


(23) 


G 


1*2 



where 

*1 = r^' 9 ' = ri (cos 0x + i sin 0x) (24) 

Splitting into real and imaginary parts and substituting for u rl and u 0l their equivalent 
expressions from equ; tions ( l) and (5) with r = r t and 0 — Oy we get 

2V r oo r i( r i “ a 2 ) 2 cos 9 1 cos(G0x + 2a) — 4V 0O rja 2 (r? - a 2 ) 2 cos0j cos40x 

+ 2V 0 Jrja 4 (r 2 - a 2 ) 2 cos 0j cos(20x - 2a) - 2V 00 rl(r* - a 4 ) sin 0 t sin(60 1 + 2a) 

+ 4V co rxa 2 (rf - a 4 ) sin 0j sin 4 0 t - 2V r 0O r 1 a 4 (rf - a 4 ) sin 0 X sin(20j — 2a) 

— 2iorJ(rJ - a 2 ) sin(60x + 2a) •+• 4 /^r}a 2 (rj — a 2 ) sin A0\ 

- 2fcorJa 4 (rJ - a 2 ) sin(20x — 2a) + 2k l r®a i sin(60x + 2a) 

+ &irja 2 (7r 2 - 11a 2 ) sin 40x + 7&xrfa 2 (r 2 — a 2 ) cos 40x 

+ 2fcjr 2 a 4 (3r 2 — 2a 2 ) sin(20x - 2a) + 6fcxr 2 a 4 (r 2 — a 2 ) cos(20x - 2a) 

+ 3(&i — mx)(r 2 — a 2 )(cos 4a — sin 4a) a° + 2mxr G a 2 cos(G0x + 2a) 

— mjrf a 2 (7r J — 3a 2 ) cos 40j — 7mxrja 2 (r 2 — a 2 ) sin 40* 

- 2mjr 2 a 4 (3r| — 4a 2 ) cos(20i — 2a) — Gmjr 2 a 4 (rJ — a 2 ) sin(20i — 2a) = 0 (25) 

ar.d 

2V 00 rf(rJ - a 2 ) 2 cos 0j sin(G0x + 2a) - 4V oc) rJa 2 (rJ - a 2 ) 2 cos0i sin 40x 

+ 2V ao ria i (r‘l — a 2 ) 2 cos 0 t sin(20i - 2a) + 2V ca r\(r\ - a 4 ) sin 0 1 cos(G0 t r 2a) 

— 4Voo rja 2 (rj — a 4 ) sin 0! cos40j 4- 2V' 0O r 1 a 4 (r 4 - a 4 ) sin 0x cos(20i - 2a) 

+ 2&orJ(rJ — a 2 )cos(G0i -f 2a) — 4fco»-Ja 2 (rJ - a 2 )cos40x 

+ 2k i) r\a i {r\ - a 2 )cos(20 t - 2a) - 2 k^a 7 cos(G0i 2a) 

- fcjrJa 2 (7rJ - 11a 2 ) cos40i — 7&xfJa 2 (rJ — a 2 ) sin 40 x 

- 2fcir 2 a 4 (3r 2 - 2a 2 )co3(20i - 2a) — GfcxrJa 4 (rJ - a 2 )sin(20x - 2a) 

— 3(/:i + mi)(r 2 - a 2 )a G (cos 4a — sin 4a) + 2m 1 r G a 2 sin(G0x 4- 2a) 

- mxr 4 a 2 (7r 2 - 3a 2 ) sin 40x - 7mxrJa 2 (rJ - a 2 ) cos40x 

— 2mxr 2 a'(3rJ - 4a 2 ) sin(20x - 2a) — Gmxr 2 a 4 (rJ - a 2 ) cos(20x — 2a) = 0 (26) 
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3.4 Non-Dimensional Parameters 


At this point some non-dimensional parameters have to be introduced if the system 
of equations (19), (21), (25), and (26) is to be solved explicitly. These are the following: 


Tf 

AO — -77— 

aV„ 

(27) 

tr _ 

A 1 — -77— 
aV'oo 

(28) 

K=-h- 

°V 00 

(29) 

aV'ca 

(30) 

ill* a. 
a 

(31) 


Substitution into the original equations yields the following system of four equations 
in five unknowns K 0 , K j, Afi, Ri, and 0i: 

(2 sin a - /f)[fl 2 + 1 - 2Ri cos^ + a)} - 2M l R l sin(0, + a) + K^R 2 - 1) = 0 (32) 
(2 sin a + K)[R\ + 1 -f 2fZj cos(0 t + a)j - 2M\Ri sin(0i + a) — Ki[R 2 - 1) = 0 (33) 

2/2j(/? 2 - l) 2 cos 0 t cos(G0j + 2a) — 4/?J(/?J — l) 2 cos 0j cos 40j 

+ 2Ri(R? - l) 2 cos 0 1 cos(20j - 2a) - 2 Ri{R* ~ 1) sin 0i sin(G0! + 2a) 

+ 4R 2 (R i — 1) sin 0j sin 40i - 2 Ri[R* - I) sin 0 1 sin (2t?i — 2a) 

- 2 K u ltf(R? - 1) sin(G0! + 2a) + - 1) sin 4 0 t 

- 2 KuR\{R\ - 1) sin (20 j - 2a) + 2 K t /l? sinf'?, + 2a) 

+ K l Ri{7Rj - U)3\n'iO l +7KiR* l {R} - l)cos i 0 X 

+ 2K l Rl(ZR 2 l -2)sin(20, - 2a) + G/f i R}{R] - l)cos(2 0, -2a) 

+ 3(/fi - Mj)(/ff - l)(cos 4a - 3in 4a) -t 2rtfj /£® cos(G0i -1- 2a) 

- M,/Z<(7rt? - 3) C03 40i - 7Af,flJ(/t? - lJsinlOx 

- 2MiR\{2R\ --4)003(20! - 2 a) - GA hR\{R\ - l)sin(20i - 2a) = 0 (3-1) 

2 R1[R\ - l) 2 co30j sin(G0i + 2a) - 4fiJ(fi 2 - l) 2 cos0 t sin 40i 

+ 2R l {R 2 l - l) 2 cos C’i sin(20! - 2 a) +2il?(/li - 1) sin 0, cos(G0 t + 2a) 


hi 



- 4R?(Rf — 1) sin 0 , cos 49 1 + 2/?,(/?f — l) sin 0 , cos(20, - 2a) 

+ 2 KoRtiRl - 1) cos(60, + 2a) - 4K 0 Rt{Rl - 1) cos 40, 

+ 2K»Rl(Rl ~ 0 cos(20, - 2a) - 2K cos (60, + 2a) 

- KiR^lR] - 11) cos 40, - 7 KiR\[R\ - 1) sin 40, 

- 2K l R\[ZR\ - 2) cos(20, - 2a) - 6 KiR\{R\ - 1) sin(20, - 2a) 

- 3(ff, + Af,)(fJ, — l)(cos 4a - sin 4a) + 2 M,i2® 3in(60, + 2a) 

- MiR^Rl - 3) sin 40, - 7A/,/^(i?? - 1) cos 40, 

- 2Af,/iJ(3i?J - 4) sin(20, - 2a) - - l)cos(20, - 2a) = 0 (35) 


3.5 Additional Condition 

We see that the equations derived so far are not enough to give us a unique solution. 
Therefore we must seek additional information in the nature of the flow. As mentioned 
in the introduction the flow is assumed inviscid which implies that the total force on the 
plate must be perpendicular to it (see Figure 4). Therefore 

tar a — — (36) 

Li 

Now we have to relate the drag and the lift with the unknowns K and A/,. 

From Figure 5 we sec that the fluid between streamlines a and b disappears into the 
sink producing drag which can be calculated \s follows: The rate of muss of fluid between 
the streamlines a and 6 (which is distance d apart) must be equal to the rate with which 
mass i3 disappearing in the sink, i.e., 

rh = —pV m db = m,p62?r. (37) 

Now the drag is given by the rate of los3 of momentum of that fluid therefore 

D = -mV*, (38) 

combining equations (37) and (38) we get 

- D 
1 2jr^6V 00 


**5 


(39) 



and from equation (30) we get 


.. \T 


but 


Afi = - 


D 

InpbaV^ 


S = be 


and from the Joukowski transformation we know that 


c = do 


also, 


C D = 


D 


yvis' 


The last four equations combine to give 


Cd — —ttA/i 


The lift can be found from the Kutta-Joukowski law 


(40) 

(41) 

(42) 

(43) 

(44) 


L = P V„ V • 6 


(45) 


and since 


and 


we get 


r = 2xk 



L 

&VZS 


Ct = vK 


(46) 

(47) 

(48) 


after employing equation (29). 

Now going back to equation (36) we sec that 

tan <* = -”• (49) 

thus wc have a system of five equations (32)-(35) and (49), which can be solved for the 
five unknowns: Ri, 0 lt K 0 , K j, and M\. This i3 done numerically (see Appendix 2) and 
the results arc shown in the next secti »n. 
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4. RESULTS 


4.1 Position of the Equilibrium Point 

As it can be seen from Figure 6 the equilibrium point is always ahead of the leading 
edge of the plate and it goes farther as the angle of attack increases reaching a maximu m 
distance of approximately 53 percent of the chord from the leading edge (measured along 
the chord line) at 45 degrees angle of attack. 

At the same time (sec Figure 7) after reaching a maximum distance above the chord 
line of approximately 16 percent of the chord (measured perpendicular to the chord line), 
it starts moving downwards crossing the chord line at a ~ 35° and getting farther under 
it at higher angles of attack. 

4.2 Vortex and Sink Strengths 

The bound vortex strength increases almost linearly with angle of attack while the 
leading edge vortex strength grows nonlinearly reaching a maximum at a = 45° and 
dropping to zero at a os 73° where the lower stagnation point reaches the trailing edge and 
moves off the plate (Figure 8). Thus the total circulation is also nonlinear with angle of 
attack and has a maximum at a — 50° (Figure 10). The sink strength on the other hand, 
starting with very small values increases almost linearly with angle of attack (Figure 9). 


4.3 Pressure Distribution on the Plate 

The pressure distribution has been calculated at the y' plane (the most convenient 
one since the plate is horizontal and lies along the £'-axisk 

From the sequence of transformations shown in Figure 2 we have 


t f u 

<=*'+- 


= ze' a + — c~ ia 


= re'°e ia + 


or 


= r [cos(0 + a) -f- t sin(0 + a)] 4- — [cos(0 + a) - »sin(0 + a)). 


^7 



I - ■> ' 


and since 


we have 


and 


f' = £' + *Y 

C ' = eos(0 + a) 

rf = sin(0 + a) 


on the surface of the plate r = a so 


t'. p = 2acos(0 + a) 

Y P = o 


(50) 

(51) 

(52) 


(53) 

(54) 


To find the velocity field around the plate, a similar procedure as in Section 2.3 yields: 

dvj _ dtu (is 
d$' ds d’f' 

dw z 1 


dz z 3 e ,a - a 2 e~ ,a 


or 


rV w 


e- K -* u O = e ’ K. -»“*)• 

and on the surface (r = u) v/c get Gnally 


1 casO 

Ur " 2 sin(0+a) U5 '* 

(55) 

1 nin0 

,tC,r r 2 sin(0 +• a) U °‘* 

(56) 

u *p = 5 u o..“c(0 + a) 

(57) 


from which 


Now from equation (5) and using equations (27)— (31) we find that 

22- = _2sin 0 — K + ^f i ^, 3 in(g-0|) + /f,(n»-l) 


1 + II l - 2Ri CO s [0 - 0 t ) 


(58) 


\ \ 
UQ 


TT/" 



Substituting equation (58) into (57) a3 well as the corresponding solution sets for each 
angle of attack we can find the pressure distribution from 


C - 1 - 

^P'P ~ 1 X/7 

w oo 


(59) 


C p , w has been plotted for four angles of attack in Figure 11 where the horizontal axis 
measures the non-dimensional parameter 


£' 

n = cosff? + a) = — ^ 
v ' 2 a 

which varies from -1 at the leading edge to +1 at the trailing edge. 

It can be seen from the plottings that both stagnation points move towards the 
trailing edge as the angle of attack increases but the lower one moves faster, being at 
about 70 percent of the chord from the leading edge at a = 60° while the upper one is 
only at the mid-chord point. 

4.4 Lift, Drag, and Pitching Moment 


At this point a comparison will be made between the results of: 


(i) the classical solution of totally attached flow over the plate with the Kutta condition 
satisfied at the trailing edge only and a singularity (infinite suction) at the lcmling 
edge (Figure 12a), 


(ii) this model with a detached vortex and a sink which make possible the satisfaction 
of the Kutta condition at botli edges removing thus the leading edge singularity 
(Figure 12b), 


(iii) the Helmholtz solution of totally separated flow over the plate (Figure 12c). 
The lift, drag and moment coefficients for the three models arc given below: 


Cr,i = 2 ~ sin a 

C12 = ~K 

7rsin2a 


Cpi = 0 


Cmi = T 3>n 
4 


C p 7 — — 5rAfj Cm 2 = ~ siu 2a 

4 


Ci 3 = 


4 + n 3in a 


Cn2 — 


2" sin 3 a 
4 + 7rsin a 


(CO) 

(Gl) 

(G2) 



From Figure 13 it can be seen that our model with partially separated Dow (as the 
leading edge vortex can be thoug.' ; of), gives the highest lift coefficient at least up to 
a = CO 9 . 

At higher angles it looks like the first model gives slightly higher lift coefficients, but 
th>3 is misleading since at high angles of attack the flow separates at some point on the 
upper surface, resulting in significantly lower lift coefficients. Thus the first model breaks 
down at high angles of attack. 


From Figure 14 it can be seen that our model gives much higher drag coefficients 
even than the third model in which the Dow is totally separated on the upper surface of the 
plate. This should not be surprising however, because the drag in our model is associated 
with a large momentum loss of all the fluid that dissappears into the sink (see Figure 5). 

As a confirmation to the above comes Figure 15 which shows the same lift to drag 
ratio for the second and third models (since L/D = tana in both cases). In other words 
since the second model yields a much higher lift than the third one, it also gives a much 
higher drag which can be thought of as induced drag (lift related drag). 


Finally, it U interesting to note that our model has the same moment coefficient and 
the same position of aerodynamic center (at the quarter chord point) a3 the first model. 
The proof is as follows: Since there is a uniform wind, the velocity at a great distance from 
the plate must tend simply to the wind velocity, and therefore if |e| is sufficiently large we 
may write 


dw 

dz 





or 


• , , ts 

w Va,z -f A In s + — H . 

z 

Now the force and the moment ca the plate can be found from the theorem of Blasius 
and it turns out (after performing the integration) that the force depends only on A while 
the moment depends only on ft. However, comparing the complex potential for the first 
two models we can sec that D is the same regardless of the presence of the vorlcx/s’mk 
combination and therefore the moment ought to he the ennic for the two models. In 
equations (GO) and (Gl) Cm is taken about the mid-point of the plate and is therefore 
positive. The same relations but with a minus sign on the right hand side arc valid for 
Cm about the leading edge. Equations (GO) and (Gl) arc plotted in Figure 1G. 
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5. CONCLUSIONS 


(i) The present inviscid, incompressible, two-dimensional model of a flat plate with a 
detached vortex close to its leading edge indicates that lift coefficients up to around 
6 are achievable. Higher values should also be possible if thickness and camber are 
added, considering an airfoil instead of a flat plate (reference 3). 

(ii) In order to satisfy the Kutta condition at both the trailing edge and the leading 
edge, the presence of the sink U necessary (see equations (32) and (33)). This b 
in agreement with physical observations of the leading edge vortex which forms 
over delta wing3 at high angles of attack where the spanwise pressure gradient due 
to sweep angle evacuates the vortex core. It is also in agreement with the model 
presented in reference 3. 

(iii) For a given lift ( C l) there are two possible solutions for the location of the equi- 
librium point and the corresponding strengths of the vortex and sink. That of the 
lower angle of attack gives a weaker sink and therefore less drag, while that of the 
higher angle of attack gives a stronger sink and the associated higher drag. The vor- 
tex strengths do not differ much for the two solutions since they arc closely related 
to the lift (which is the same for the two solutions). This result b also in agreement 
with reference 3. 

(iv) The upper limit found in the Cl versus a curve suggests that if a stronger vor- 
tex would exist at the' equilibrium point, unrealistic supcrcirculation would occur 
resulting in the streamlines going entirely around the system. 

(v) A limitation of the present model appears at n ~ 73°, above which the lower stag- 
nation point moves ofT the plate making thus impossible the Qow pattern depicted 
in Figure 5, on which thb model is based. 

(vi) A comparison between Figures 8 and lc shows very good agreement (at least qualita- 
tively) between the results derived here and those found experimentally in reference 

1 . 

(vii) Finally, the presence of the vortex and the sink docs not affect the position of the 
aerodynamic center which remains at the quarter chord point. 




»*r 
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Schematic of the Dow over a slender sharp-edged wing. 




Figure 1. 







c. Joukowski plane (f'-plane). 


d. Pinal plate plane (f-plane). 


Steps in mapping sequence from original circle plane to final flat airfoil 
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Figure 7. Position of the equilibrium point (perpendicular to the plate). 
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Figure 8 . 


Leading edge and bound vortex strength versus angle of attack. 
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APPENDIX 1 


A1.1 

Consider a vortex of strength ki at the point *i outside of a cylinder |z| = a. Then, 
the complex potential is . 

w = /(*) + 1 ( 7 ) 

= iki ln(z — ri) - »&i In z^j 

according to the circle theorem (reference 2) w can also be written as 

w = iki In z + iki ln(r — z t ) — iki In + const. 

which shows three vortices: one at the point «i; one at the point ^ which is the inverse 
square point of zj with respect to the cylinder; and one at the origin (center of the cylinder). 

A1.2 

Consider a sink of strength m j at the point r t outside of a cylinder |z| = a. Then, 
the complex potential is 

« = /(*)+/( 7) 

w = mj ln(z — zj) + mi In £ 1 ^ 

or 

w = — mi lnz +m t ln(z - z{) + m t In +• const, 

equivalent to: a sink at the point * 1 , a sink at the point and a source at the origin. 
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APPENDIX 2 


A four-dimensional Newton-Raphson algorithm has been used to solve the system 
of equations (32), (33), (34), and (35) for R j, 0 1 , /f 0 , K i after elimination of iVfj from 
equation (49). 



n n o 


PROGRAM NR4D 

REAL K0.K1.K.J2.J3.J4.J5,J6.J7, II, 12. 13, 14, IS, 16,17.18,19 


DATA 


ALFA=? 

ALF=ALFA* 3 . 14159/180 . 0 
30 WRITE (5.110) 

110 FORMAT ( * GIVE UNDERRELAXATION FACTOR*) 

READ (5,120) C 
120 FORMAT (F) 

C T=TnETA, R=R1 

R=? 

T=? 

K0=? 

Kl=? 

10 WRITE (5,100) R,T,K0,K1 

100 FORMAT (4E15.6) 

C 

T1=6*T*2‘ALF 

T2=2‘T-2*ALF 

T3=T+ALF 

K=K0*K1 

C 

A1=2*R* *4* (R‘R-1) * (9*R*R-S) *COS(T) ‘COS(Tl) 

A2=-4*R*R* (R‘R-1) * (7*R‘R-3) *COS (T) ‘COS (4‘T) 

A3=2* (R‘R-1) * (5*R*R-1) ‘CCS (T) ‘COS (T2) 

A4=-2*R“4* (9*R* *4-5) »SIN(T) *SIN(T1) 

A5=4*R*R* (7*R* *4-3) *SIN(T) *SIN(4‘T)-2* (5‘R“4-1) *SIN(T) ‘SIN (72) 
A6=-4*K0*R“5* (4»R*R-3) ‘SIN(Tl) +8*K0*R‘ *3* (3‘R*R-2) »SIN(4*T) 
A7=-4‘K0*R» (2‘R‘R-l) »SIN(T2) + 12*K1‘R“5‘SIN (Tl) 

A8=2*K1*R“3* (21*R*R-22) »SIN(4*T) +14*K1*R“3* (3‘R»R-2) *C0S(4»T) 
A9=8*K1*R* (3*R‘R-1) *SIN(T2) *12*K1*R* (2*R*R-1) ‘COS (T2) 

B1=6*R* (K1+K‘TAN (ALF) ) * (COS (4*ALF) -SIN(4*ALF) ) 
B2=-12‘R**5*K‘TAN(ALF) *COS (Tl) 

B3=3*K*R“3‘ (14*R*R-4) ‘TAN (ALF) »COS(4*T) 

B4=14‘K*R“3* (3*R*R-2) ‘TANfALF) *SIN(T) 

B5=K*R*TAN (ALF) * ( 8 » (3*R»R-2) *COS(T2)+12* (2‘R‘R-l) *SIN(T2)) 

C 

All = Al + A2 + A3+A4+A5 + A6*A7+A8+A9+Bl + B2-*-B3+B4+B5 

C 

B6=-2*R“5‘ (R‘R-1) * *2* (SIN(T) ‘COS (Tl) -*-6*C0S (T) *SIN(T1)) 
B7=4‘R“3* (R‘R-1) “2* (SIN(T) ‘COS (4*T) +4*C0S (T) *SIN(4*T) ) 
B3=-2‘R» (R*R-1) “2* (SIN(T) *C0S(T2) +2‘COS(T) *SIN(T2) ) 

B9=-2*R“5* (R* *4-1) * (COS (T) *SIN (Tl) -* 6 ‘SIN (T) ‘COS (Tl) ) 

C1=4*R**3* (R* * 4 - 1 ) • (COS (T) *SIN(4*T) +4‘SIN(T) *COS(4*T)) 

C2=-2‘R* (R“4-l) * (COS (T) *SIN(T2) +2*SIN(T) *COS(T2) ) 

C3=-12*K0*R* * 6 * (R*R-1) “COS (Tl) ♦16*K0*R“4‘ (R‘R-1) *C0S(4*T) 
C4=-4*K0*R*R* (R‘R-1) ‘COS (T2) +12*K1*R“6‘C0S (Tl) 

C5=4‘K1‘R“4* (7*R*R-11) ‘CCS (4*T) -2Q*K1*R“4* (R‘R-1) *SIN(4‘T) 
C6=4*K1»R*R* (3*R*R-2) ‘COS (T2) -12*K1«R*R* (R‘R-1) *SIN(T2) 
C7=4*K»R“4»TAN(ALF) * (3*R‘R‘SIN(T2) -(7*R‘R-3) *SIN(4*T)) 
C8=28‘K*R“4* (R“2-l) ‘TAN (ALF) *C03(4»T) 

C9=-4*K*R*R*TAN(ALF) • ( (3*R*R-4) ‘SIN (T2) -3* (R‘R- 1 ) ‘COS (T2) ) 


A12 = E5 + B7+5S + B9+C1-»C2*C3+C4+C5+C6+C7+C3+C9 

D1=-2*R“6* (R‘R-1) *SINfTl)+4*R“4* (R‘R-1) *SIN(4*T) 
D2=-2*R*R* (R‘R-1) *SIN(T ’) 



ORIGINAL PAGE (S 
OF POOR QUALITY 


D3=3* (R*R-1) *TAN (ALE) * (COS (4*ALF) -SIN(4‘ALF) ) 

D4=-2*R**6*TAN (ALT) *COS (Tl) +R* *4* (7*R*R-3) *TAN (ALF) *CCS (4*T) 
D5=7*R**4* (R*R-1) *TAN (ALF) *SIN(T) 

D6=2*R*R* (3*R*R-4) *TAN (ALF) *COS (T2) 

D7=6*R*R*(R*R-1) •TAN(ALF) *SIN(T2) 

A13=D1+D2*D3+D4*D5*DS*D7 

D8=2*R* *6*SIN (Tl) +R**4* (7«R*R-U) *SIN (4*T) 

D9=7*R**4* (R*R-1) *C0S(4*T) *2*R*R* (3*R*R-2) *SIN(T2) 

E1=6*R*R* (R*R-1) *COS (T2) 

E2=3‘ (1+TAN (ALF) ) * (R*R-1) * (COS (4*ALF) -SIN (4*ALF) ) 
E3=-2*R**6*TAN(ALF) *COS (Ti)*R»*4* (7»R*R-3) *TAN (ALF) »COS (4*T) 
E4=7*R**4* (R*R-1) *TAN(ALF) *SIN(T) 

E5=2*R*R*TAN (ALF) * ( (3*R*R-4) *COS (T2) +3* (R*R-1) *SIN(T2) ) 
A14=D8+D9+El*E2+E3*E4+E5 

E6=2*R* *4* (R*R-1) * (9*R*R-5) *COS(T) *SIN(T1) 

E7=-4*R*R* (R‘R-1) * (7*R»R-3) *COS(T) *SIN(4*T) 

E8=2* (R*R-1) * (5*R*R-1) *CCS(T) *SIN(T2) 

E9=2*R* *4* (9*R“4-5) *SIN(T) *COS(Tl) 

G1=-4*R*R* (7*R»*4-3) »SIN(T) *COS(4*T)+2* (5*R**4-1) *SIN(T) *COS (T2) 
G2=4»K0*R*»5* (4*R*R-3) *CCS (Tl) -8*KD*R**3* (3*R*R-2) *COS (4*T) 
C3=4*K0*R* (2*R*R-1) *CCS (T2) -12*K1*R* *5*COS (Tl) 

G4=-2*K1*R* • 3* ( (21*R*R-22) »COS (4*T) * 7 * (3*R*R-2) *SIN(4*T) ) 
G5=-4*K1*R* (2* (3*R*R-1) *COS(T2) +3* (2*R*R-1) *SIN(T2)) 

G6=-6* (K1-K*TAN (ALF) ) *R» (COS (4*ALF) -SIN (4*ALF) ) 
G7=-12*K‘R«*5*TAN(ALF) *SIN(T1) 

G8=K*R* *3*TAN (ALF) * (3* (14*R*R-4) »SIN(4*T) *14* (3*R*R-2) 'C0S(4*T)) 
G9=4*K*R*TAN(ALF) * (2* (3*R*R-2) *SIN(T2)*3» (2*R*R-1) *C0S(T2)) 

C 

A21 = E6*E7*E8+E9+G1*C2*G3+G4+GS*G6*G7+C8+G9 

C 

I1=2*R**5* (R*R-1) **2* (-SIN (T) *SIN (Tl) *6*COS (T) *COS (Tl) ) 
I2=-4*R**3* (R*R-1) **2* (-SIN (T) *SIN (4*T) +4*CCS (T) *C0S(4*T)) 
I3=2*R* (R»R-1) * *2* (-SIN (T) *SIN(T2) *2»COS(T) *COS(T2)) 

14=2 *R* *5* (R* *4-1) * (COS (T) *COS (Tl) -6* SIN (T) *SIN (Tl) ) 

IS=-4*R**3* (R**4-l) * (CCS (T) »COS (4*T) -4*SIN(T) *SIN(4*T)) 

I6=2*R* (R* *4-1) * (COS (T) *COS (T2) -2* SIN (T) *SIN (T2) ) 

I7=-4*K0*R*R* (R*R-1) * (3*R**4‘SIN(T1) +SIN(T2) -4»R*R*SIN(4*T) ) 
I8=4*K1*R* *4* (3*R*R*3IN(T1)*(7*R*R-11) *SIN(4*T)) 

I9=-4*K1*R*R* (7*R*R* (R*R-1) *COS (4*T) - (3*R*R-2) *SIN(T2)) 
P1=-12*K1*R*R* (R*R-1) *COS (T2) 

P2=-4*K*R* *4*TAN (ALF) * (3*R*R*COS (Tl) - (7*R*R-3) *COS (4‘T) ) 
P3=-28*!(*R* *4* (R* *2-1) * TAN (ALF) *SIN(4*T) 

P4=4*K»R*R*TAN(ALF) * ( (3*R*R-4) *COS(T2)-3* (R*R-1) *SIN(T2) ) 


A22 = Il+I2*I3*I4*I5+I6*I7*IG*I9*Pl+?2+P3*P4 

P5=2*R*R* (R*R-1) * (R« *4«COS (Tl) -2*R*R*COS (4*T) +COS (T2) ) 
P6=3* (R*R-1) "TAN (ALF) * (COS (4*ALF) -SIN (4*ALF) ) 
P7=-2*R**6*TAN(ALF) *SIN(T1) 

P8=R* *4*TAN (ALF) * ( (7*R*R-3) *SIN(4*T) *7* (R*R-1) »C0S(4*T) ) 
P9=2*R*R*TAN (ALF) * ((3*R*R-4) *SIN(T2) +3* (R*R-1) *COS (T2) ) 


A23=P5*P6+P7*P8+P9 


C 

Ql=-R* *4* (2*R°R*COS (Tl) * (7*R*R-11) *COS (4*T) *7* (R*R-1) *SIN (4*T) ) 
Q2=-2*R»R» ( (3*R*R-2) *CC3 (72) *3* (R*R-1) *SIN(T2) ) 




c- 

c 


c 


c 


c 

c 


c- 


Q3=-3* (1-TAN (ALF)) * (R*R-1) * (COS (4*ALF) -SIN (4*ALF) ) 
Q4=-2»R»*6*TAN(ALF) *SIN(T1) 

Q5=R**4*TAN(ALF) * ( (7*R‘R-3) *SIN(4*T) +7* (R»R-1) *C0S(4*T) ) 
Q6=2*R*R*TAN (ALF) * ( (3*R‘R-4) *SIN(T2) +3* (R»R-1) *COS(T2) ) 


A24=Q1 + Q2 + Q3+Q4 + Q5*Q6 


A31=2* (2*SIN (ALF) -K) * (R-COS (T3) ) +2*K1‘R+2*K*TAN(ALF) *SIN(T3) 
A32=2*R* (2*SIN (ALF) -K) »SIN(T3) *2*K*R*TAN(ALF) *COS (T3) 

A33=- (l-*R*R-2*R*COS (T3) ) +2*R*TAN (ALF) *SIN(T3) 

A34=- (l + R*R-2*R*COS (T3) ) +R*R-1+2‘R*TAN (ALF) »SIN(T3) 


A41=2* (2*SIN (ALF) +K) * (R+COS (T3) ) +2*K‘TAN (ALF) *SIN (T3) -2*K1‘R 
A42=2*K*R*TAN (ALF) ‘COS (T3) -2»R* (2*SIN(ALF) +K) *SIN(T3) 
A43=l+R* *2 + 2*R*CCS (T3) -*-2*R*TAN (ALF) *SIN(T3) 

A44=l+R* * 2 + 2*R*COS (T3) + 2*R*TAN (ALF) *SIN (T3) - (R*R-1) 


FA=2*R* (R*R-1) **2*COS(T) * (R“4*COS (Tl) -2*R*R*COS (4‘T) +COS(T2) ) 
FS=-2*R* (R“4-l) *SIN(T) * (R“4*SIN(T1) -2»R»R*SIN(4‘T) +SIN(T2) ) 
FC=-2*KO*R*R* (R*R-1) * (R“4‘SIN(T1) -2‘R*R*SIN(4‘T) +SIN(T2) ) 
FD=2*K1*R* *6* SIN (Tl) 

FE=K1*R* *4* ( (7*R*R-11) ‘SIN(4*T)*7* (R»R-1) »C05(4»T)) 
FG=2*K1*R*R* ( (3*R*R-2) *SIN(T2) *3* (R*R-1) *COS(T2)) 

FH=3* (K1*K*TAN (ALF) ) * (R‘R-1) • (COS (4*ALF) -SIN (4*ALF) ) 
FI=-2*K*R“6*7AN(ALF) »COS(Tl) 

FJ=K*R**4*TAH(ALF) * ((7*R*R-3) *COS(4*T) ♦7* (R*R-1) *SIN(4*T)) 
FK=2*K*R*R*TAN (ALF) * ( (3*R‘R-4) ‘COS (T2) +3* (R*R-1) *SIN(T2) ) 

FI = FA+FB+FC+FD*FE+FG+FH-»FI*FJ‘-FK 


FL=2*R* (R*R-1) “2‘CCS (7) * (R“4*SIN(71) -2*R fi R*SIN (4*1) +SIN(T2) ) 
FM=2‘R* (R* *4-1) ‘SIN (T) * (R“4*C0S (Tl) -2*R‘R‘COS (4*T) -COS (T2) ) 
FN=2‘KO‘R*R* (R*R-1) * (R“4*COS (Tl) -2*R*R*COS (4*7) +COS (T2) ) 
F0=-2*K1*R“6 4 C03 (Tl) 

FP=-K1*R* *4* ( (7*R*R-11) *COS (4*T) -*-7* (R*R~1) *SIN (4*T) ) 
FQ=-2*K1*R*R‘ ((3‘R‘R-2) *COS(T2)+3‘ (R‘R-1) *SIN (T2) ) 

FR=-3* (iCl-K*TAN (ALF) ) * (R»R-1) * (CCS (4*ALF) -SIN (4*ALF) ) 
FS=-2‘K*R**6‘TAN(ALF) *SIN(T1) 

FT=K*R“4*TAN (ALF) * ( (7*R*R-3) *SIH (4*T) +7‘ (R*R-1) *COS (4*T) ) 
FU=2*K*R*R*TAN (ALF) • ((3*R*R-4) ‘SIN (T2) +3* (R*R-1) ‘CCS (T2) ) 


F2 = FL + EM+FN+FO-*-F?+FQ+FR + FS+FT*FU 

C 

FV=(2*SIN(ALF)-K) * (1*R*R-2*R*CCS (T3) ) 

FW=K1* (R* *2-1) +2*K*R*TAN (ALF) *SIN(T3) 

F3=FV*FW 

FX= (2*SIN (ALF) ■‘•K) * (l+R e R+2 *R‘COS (T3) ) 

FY= 2 * IC * R * TAN (ALF ) *SIM(T3)-K1* (R*R-1) 

F4=FX*FY 

S1=A11*A22* (A33*A44-A43*A34) -*-All*A23* ( A3 2 » A44-A4 2 * A3 4 ) 
S2=A11*A24* (A32*A43-A42‘A33) +A12"A21* (A33*A44-A43*A34) 
S3=A12*A23* (A31‘A44-A34*A41) +A12*A24* (A31*A43-A41*A33) 
S4=A13*A21* ( A.32 * A44-A34 * A42 ) +A13»A22* (A31*A44-A34«A41) 
S5=A13*A24* (A31 * A42-A41 *A32 ) +A14*A21* (A32*A43-A33*A42) 

C 6=A14*A22* (A31»A43-A33*A41) +A14»A23“ (A31*A42-A41*A32) 

DET = S1+S2+S3+S4+S5+S6 


J2=F1‘A22* (A33*A44-A43‘A34) +F1*A23* (A32*A44-A42*A34) 
J3=F1*A2 4* (A32*A43-A42*A33) +A12‘F2* (A3 3 * A44-A4 3 * A34) 
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J4=A12*A23* (F3*A44-A34*F4) +A12*A24* (F3*A43-F4*A33) 
J5=A13*F2* (A32*A44-A34*A42) +A13*A22» (F3*A44-A34*F4) 
J6=A13*A24* (F3*A42-F4*A32) ■»A14*F2* (A32*A43-A33*A42) 
J7=A14*A22* (F3*A43-A33*F4) *A14*A23* (F3*A42-F4*A32) 

HI = (J2-J3-J4>J5+J6»J7) / DET 

U1=A11*F2* (A33*A44-A43*A34) +Aii»A23» (F3*A44-F4*A34) 
U2=A11*A24* (F3*A43-F4*A33) +Fl r A21* (A33*A44-A43* A34) 
U3=F1*A23* (A3i*A44-A34*A41) +F1*A24* (A31»A43-A41»A33) 
U4=A13*A21* (F3*A44-A34*F4) ♦A13«F2* (A31*A44-A34»A41) 
U5=A13*A24* (F4*A31-F3*A41) +A14*A21* (F3*A43-F4*A42) 

' U6=A14*F2* (A31*A43-A33*A41) -»A14*A23* (F4*A31-F3*A41) 

H2 = (U1+U2-*U3+U4*U5*U6) / DET 

V1=A11*A22* (F3*A44-F4*A34) ♦A11*F2* (A32*A44-A42*A34) 
V2=A11*A24* (A32*F4-A42*F3) t-AI2 # A21* (F3*A44-F4*A34) 
V3=A12*F2* (A31*A44-A34*A41) +A12»A24» (A31*F4-A41*F3) 

V4=F1 *A21* (A32*A44-A34*A42)«-F1*A22* (A31 «A44-A34»A41) 
V5=F1*A24* (A31*A42-A41«A32) t-A14*A21* (A32*F4-F3*A42) 
V6=AI4*A22* (A31*F4-F3*A41) *A14*F2 e (A31*A42-A41«A32) 

H3= (V1+V2+V3+V4+V5+V6) / DET 

C 

W1=A11*A22* (A33*F4-A43*F3) «-All*A23» (A32*F4-A42*F3) 
W2=A1*F2* (A32*A43-A42*A33) +A12*A21* (A33*F4-A43*c3) 
W3=A12*A23* (A31*F4-A41»F3) +A12*F2* (A31*A43-A41*A33) 
W4=A13*A21* (A32*F4-A42*F3) +A13*A22* (A31*F4-A41*F3) 
W5=A13*F2* (A31*A42-A41*A32) +F1*A21* (A32*A43-A33*A42) 
W6=F1*A22* (A31*A43-A33*A41) *F1*A23* (A31*A42-A41*A32) 

K4 = (W1+W2+W3+W4+W5+W6)/ DET 

C 

R =R-C»H1 
T =T-C*H2 
K0=K0-C*H3 
K1=K1-C*H4 

IF (ABS (HI) .LT. 1 .E-3 .AND. ASS (H2) .LT. 1 .E-3) GO TO 20 
IF(A3S(H3) ,LT.1.E-3.AND.ABS(H4) .LT.l.E-3) CO TO 20 
CO TO 10 

20 WRITE (5, ICO) R.T.K0.K1 

WRITE (5,130) 

130 FORMAT ( * IF YO’J WANT ANOTHER RELAXATION FACTOR. TYPE 1*) 
READ (5,140) IFF 

140 FORMAT (I) 

IF (IFF.EQ.l) CO TO 30 

STOP 

END 
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